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PREFACE 


This report describes the work performed from December 5, 1978, 
through July 5, 1980, by the McDonnell Douglas Research Laboratories 
(MDRL) under the National Aeronautics and Space Administration (NASA) 
Contract NAS8-33107. The ultimate objectives of the study are to 
quantitatively establish the characteristics of Hg^ ^Cd^Te as grown on 
Earth as a basis for subsequent evaluation of the material processed In 
space and to develop the experimental and theoretical analytical 
methods required for such evaluation. The technical report was 
submitted to NASA In July 1980. 

The work was performed In the Solid State Sciences Department, 
managed by Dr. C. R. Whltsett. The Principal Investigator was 
Dr. S. L. Lehoczky; Co-Investlgators were Dr. F. R. Szofran and 
Dr. B. G. Martin. 
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1. INTRODUCTION 



1 . 1 Scope of Study 

The overall study includes the entire range, 0 < x < 1, of Hgi-x^dj^Te 
alloy compositions. Crystals were prepared by the Bridgman method with a wide 
range of crystal growth rates and temperature gradients adequate to prevent 
constitutional supercooling under diffusion-limited, steady-state, growth 
conditions. The longitudinal compositional gradients for different growth 
conditions and alloy compositions were calculated and compared with 
experimental data to develop a quantitative model of the crystal growth 
kinetics for the Hgj_jjCdjjTe alloys, and measurements were performed to 
ascertain the effect of growth conditions on radial compositional gradients. 
The pseudobinary HgTe-CdTe constitutional phase diagram was determined by 
precision differential-thermal-analysis measurements and used to calculate the 
segregation coefficient of Cd as a function of x and interface uemperature. 
Experiments will be conducted and theoretical analysis will be performed 
during the next phase of the program to determine the ternary phase equilibria 
in selected regions of the Hg-Cd-Te constitutional phase diagram, and electron 
and hole mobilities as functions of temperature will be analyzed to establish 
charge-carrier scattering probabilities. Computer algorithms specific to 
Hgi-xCdxTe were developed for calculations of the charge-carrier 
concentrations, charge-carrier mobilities. Hall coefficient, optical 
absorptance, a 'd Fermi energy as functions of x, temperature, ionized donor 
and acceptor concentrations, and neutral-defect cor.t ent rat ions . 

1.2 Summary of Progress 

A series of differential-thermal-analysis (DTA) measurements was 
performed for Hg^_jjCdjjTe alloy compositions with x ■ 0, 0.1, 0.2, 0.3, 0.4, 
0.6, 0.7, 0.8, 0.9, and 1.0. The llquldus and solidus temperatures deduced 
from the DTA data were used to establish the pseudobinary HgTe-CdTe 
constitutional phase diagram. The segregation coefficient of Cd was 
determined as a function of x and Interface teiq>erature. 

Iterative phase-equilibria calculations based on the regular associated 
solution (R.A.S.) theory^ were performed, and a set of R.A.S. parameters was 
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obtained by simultaneously least-squares-f itting the binary Hg-Te and Cd-Te^~^ 
and pseudobinary llgTe-CdTe phase diagrams. The R.A.S. parameters were used to 
calculate the activities of Kg, Cd, and Tc 2 and their partial pressures over 
pseudobinary melts. The calculated mercury partial pressures are in 
reasonable agreement with Stelnlnger's experimental results.^ The validity of 
the R.A.S. parameters and thus the calculated partial pressures are 
sensitively dependent upon necessary assumptions in the calculations; for 
example, negligible Te-Te association in the melt was assumed. 

Hgo^8Cdo.2^® alloy crystals were grown by the Bridgman method at constant 
furnace translation rates of 0.0685, 0.310, 0.597, 1.12, and 5.62 um/s. The 
temperature gradients at the solid-liquid interface, except for the 5.62 pm/s 
growth rate, were adequate to prevent constitutional supercooling in the 
melts. For three of the ingots, the longitudinal compositional profiles were 
determined by precision density measurements and were compared with calculated 

O 

profiles" for various assumed values of D, the liquid HgTe-CdTe Interdiffusion 
coefficient. The calculations included the variations in the interface 
translation rates and Interface segregation coefficients during the growth 
process. These variations were calculate, from the temperature profiles 
measured prior to each growth run, the measured compositional profiles, and 
phase equilibria data. The best fits to the data for the alloys with x = 0.2 
yielded 5.5 x 10“^ cm^/s for D and 1.5 x 10^ K*s/cm^ for G/R, the interface- 
temperature-gradient to growth-rate ratio required to prevent constitutional 
supercooling. Radial compositional variations were measured on thin slices 
from several Ingots by Infrared (IR) transmission-edge mapping. The radial 
compositional profiles deduced from the cut-on wavelengths implied concave 
solid/liquid Interfaces for the entire lengths of the crystals. 

Theoretical models and computer programs specific to Hgj^_j^CdjjTe were 
developed for calculations of charge-carrier concentrations. Hall c' ?f f iclent , 
Fermi energy, and conduction-electron mobility as functions of x, temperat re, 
and ionized-defect and neutral-defect concentrations. As in previous work on 
the HgCdSe alloy system^* the Kane three-band model^^ was used to describe 
the band structure of the HgCdTe alloys, and the best available band 
parameters were comr iled from the literature. The temperature dependence of 
the electron concentration was calculated for various net donor concentrations 
from 4.2 K to 300 K, and the calculated results agreed well with available 
data . 



The mobility calculations Included the following Instrlnsic scattering 
mechanisms: longltudlnal-optlcal phonon, longitudinal and transverse- 

acoustical phonon, heavy-hole, and alloy disorder potential. The extrinsic 
scattering mechanisms Included charged and neutral point-defects. A 
comparison of calculated results with available experimental data Indicated 
that longitudinal optical-phonon and charged and neutral defect scattering are 
the dominant mobility-limiting mechanisms. 



PSEUDOBINARY HgTe-CdTe PHASE DIAGRAII 


2.1 Alloy Preparation 

The HgTe-CdTe alloys used for the phase diagram determinations were 
prepared by reacting the constituent elements in sealed, f used-silica, 5-nnn 
i.d. X 10-mm o.d. ampules. The ampules were cleaned and etched in aqueous HF 
and annealed at 1150°C in vacuum to remove residual contaminants. Square 
cross-section bars of 99.9999% pure Cd and Te were cut from large Ingots. The 
Cd bars were etched in concentrated HNO^ and repe. cedly rinsed in methanol. 

The Te bars were etched in Br2 and repeatedly rinsed in methanol. The ampules 
were loaded in a vertical position first with 99.99999% pure Hg, then Te, and 
finally Cd. This procedure prevented contact between Hg and Cd in the 
presence of air. The ampules were evacuated and backfilled with He several 
times before final evacuation and sealing. Using the above procedure, little 
or no wetting occurred between the ampule and the alloy during reaction of the 
elements. In many cases, the reacted alloy would slide Inside the ampule. 

The various alloy preparations used for the phase diagram determinations 
are summarized in Table 1. As implied by the table, a substantial number of 
prepared ampules were discarded for var'ous reasons. Following a series of 
initial experimental runs, a concentra'.ed effort was undertaken to increase 
the fill factor (reduce the f ree-volumc) of the ampules. For this purpose, a 
series of alloys was cast, the ampules were reopened, solid quartz rods were 
placed in the ampules, and the ampules were resealed. For the ampules 
containing alloy compositions of x = 0. 1 and 0.2, the rod diameters were 
4 mm. For alloy ompositions with x-valuts of 0.3 to 0.9, 5-mm diameter rods, 
etched in concentrated hydrofluoric acid to achieve close fits, were used. 

With the exceptions of the 0.1 and 0.2 alloys, this sealing procedure 
repeatedly resulted in ampule ruptures at the measurement temperatures, and 
the cause of these failures is unclear. 

In subsequent preparations of alloys in the hlgh-x range, an increase of 
the alloy volume in the ampule was used to Increase the ampule fill-factor. 

The alloy constituents for the ampules that served as samples for the majority 
of the measurements described in this report are listed in T'fble 2. 
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TABLE 1. ALLOY SAMPLES PREPARED FOR THE PHASE-DIAGRAM 
DETERMINATIONS. 



Sample 

no. 

Com position 
(x) 

Quantity 

(moles) 

1 

Comments 

1 

0.1 

0.02 

A, B 

2 

0.2 

0.02 

A, B 

3 

0.2 

0.06 

A, C 

4 

0.1 

0.06 

A, C 

5 

0.4 

0.06 

A, B, C 

6 

0.3 

0.06 

A,C 

7 

0.5 

0.06 

A 

8 

0.6 

0.06 

A 

9 

0.7 

0.06 

A, B 

10 

0.8 

0.06 

A. B 

11 

0.9 

0.06 

B 

12 

0.8 

0.06 

B 

13 

0.2 

0.06 

C 

14 

0 

0.06 

A, C 

15 

0 

0.06 

D 

16 

1.0 

0.03 

0 

17 

0.7 

0.06 

B 

18 

0.5 

0.06 

C, F 

19 

0.7 

0.06 

C, F 

20 

0.8 

0.06 

C,F 

21 

0.1 

0.06 

D 

22 

0.2 

0.06 

0 

23 

0.3 

0.06 

C, F 

24 

0.4 

0.06 

D 

25 

0.9 

0.06 

C, F 

26 

0.8 

0.06 

C 

27 

0.9 

0.06 

C 

28 

0.9 

0.12 

D 

29 

0.8 

0 12 

D 

30 

0.7 

0.12 

C 

31 

0.6 

0.12 

C 

32 

0.6 

0.12 

C 

33 

0.7 

0.12 

C 

34 

0.5 

0.12 

C 

35 

0.3 

0.06 

D 

36 

0.5 

0.06 

C 

37 

0.6 

0.06 

C 

38 

0.5 

0.06 

C 

39 

0.6 

0.06 

D 

40 

0.7 

0.06 

0 


* A - Used for preliminary experiments D ■ Used for measurements reported here 

B ■ Some sticking of alloy to ampule E - 0.4 mm diam rod sealed inside after reacting alloy 
C - Ampule burst F ■ As E, except that an HF-etched'0.5-mm rod vyas used 
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TABLE 2. ELEMENTAL CONSTITUENTS OF ALLOYS USED FOR THE PHASE 
DIAGRAM DETERMINATION. MASS MEASUREMENT ACCURACY 
IS :t0.56 ms- 


Composition 

M 

Hg 

(g) 

Cd 

(g) 

Te 

(g) 

Stoichiometric 

Hg 

(g) 

Excess 

Hg 

(mg) 

0.1 

10.8349 

0.6732 

7.6555 

10.8345 

0.4 

0.2 

9.6299 

1.3501 

7.6570 

9.6289 

1.0 

0.3 

9.5487 

2.2950 

8.6748 

9.5428 

5.9 

0.4 

7.2201 

2.6992 

7.6555 

7.2190 

1.1 

0.6 

4.5719 

3.8343 

7.2582 

4.5687 

3.2 

0.7 

3.6452 

4.7582 

7.7183 

3.6433 

1.9 

0.8 

4.7977 

10.7427 

15.2430 

4.7939 

3.8 

0.9 

2.40S3 

12.1474 

15.3200 

2.4083 

1.0 


6P0M7S0-3 


2.2 Experimental Method for Differential Thermal Analyses 

Differential thermal analysis (DTA) measurements were used for the 
determination of the phase equilibria temperatures. The experimental 
arrangement of the ampules in the DTA furnace is shown in Figure 1. 

An ampule containing the alloy sample is mounted end-to-end with a second 
ampule coaxially inside the 460-mm long isothermal furnace liner (Na heat 
pipe) and the 610-mm long tube furnace. The second ampule contains antimony, 
which in addition to serving as the reference heat capacity provides a 
thermocouple calibration point for each experimental run. The ampules are 
held on axis by a quartz centering fixture (not shown in the figure). The 
annular volume between the ampules and a 25-mm o.d. quartz tube is filled with 
powdered diatomaceous earth. The 25-mm tube is inside an Inconel tube, which 
prevents damage to the heat pipe and furnace in case of ampule failure. 

Both the alloy and Sb-re£erence ampules are fitted with a Pt-foil band 
~5-mm wide and positioned near the middle of the ampule axially. The two 
junctions of a type-K differential thermocouple are secured to the bands. A 
separate cype-K thermocouple is secured to the band on the Sb-reference 
ampule. Each thermocouple lead extends to a copper junction maintained in a 
trlple-point-of-water cell. A block diagram of the furnace control and data 
acquisition instrumentation is shown in Figure 2. The temperature of the 
Sb-reference sample was monitored on the x-axls, and the differential 
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Figure 1. Experimental arrangement for differential thermal analytis meaturements. 

thermocouple output on the y-axis of an x-y recorder. Each graph— paper used 
for data recording was calibrated by drawing vertical lines at 1— mV intervals 
by driving the x-axis with a potentiometer (Leeds and Northrop model 8691). 
This potentiometer was in turn periodically checked against another 
potentiometer (Leeds and Nortnrop K— 4) and verified to be accurate to within 
±20 uV. The output of a type-K thermocouple is approximately 40 pV/deg. 
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The temperature programmer allows the furnace temperature to be changed 
at a uniform rate and turns the x-y recorder servo and pen-lift on and off. 
The furnace protection circuit, a meter relay which drives a power relay, 
provides an additional margin of safety against excessive furnace 
temperatures . 

A thermocouple calibration facility, similar to that described in NBS 
Circular 590^^, using graphite crucibles and NBS reference materials was 
constructed and used to calibrate thermocouple wire lots used for the 
measurements. The melting-temperature reference materials were Zn, Al, and 
Cu. The melting temperature of t' ^ Sb-reference sample served as an internal 
calibration temperature, and the Sb melting point was checked at least once 
prior to each experimental run. 

2.3 Differential Thermal Analysis Results 

Generally, for DTA measurements, the output signal from the differential 
sensors varies directly with the heating/cooling rates used, i.e., high 
heating/cooling rates yield larger signals and thus better signal to noise 
ratios. Measurements made at too high heating/cooling rates, however, can 
significantly alter the details of the DTA curves, thereby yielding phase 
transformation temperatures that do not reflect true equilibrium values. To 
assure that heating/cooling rates used in this study were within the allowable 
range, HgTe calibration sample test runs were made at a series of heating 
rates (0.5, 1, 2, 3, 4, and 5°C/min) . The thermal arrest curves for the 
various runs are shown in Figure 3. The melting temperatures deduced from the 
heating curves are within ±1.7°C of 669. 5°C for the different heating rates, 
and there is no evidence for systematic variation of the onset of thermal- 
arrest with heating rate. Because the maximum cooling rate of the furnace at 
or below 670°C was about 3°C/mln, all alloy heating and cooling curves were 
recorded at a rate of 2°C/min in order to enhance the differential 
thermocouple outputs while allowing Identical heating and cooling rates. 

Alloy Inhomogeneities broaden the range of temperatures over which 
transformations otcur and thus introduce significant uncertainties into the 
phase equilibria temperatures as deduced from DTA. For this reason, the 
samples were always rapidly quenched following casting and uniformly cooled 
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Chromei-alumel thermocouple voltage (mV) 


24 25 26 27 28 29 30 31 32 



Figure 3. Rate dependence of thermal arrest for HgTe. 


after each DTA run to reduce the possibility of preferential alloy segregation 
along the a ipule axis* To eliminate radial alloy Inhomogenelties, the samples 
were annealed at temperatures, T^, just below their solidus values • To 
awe*" tain the minimum annealing times required, samples with x » 0*1, 0,2, 0.8 
ana 0.9 were annealed for successively longer annealing times from 1 to 60 h, 
and a DTA run wa- made following each anneal. An increase in the annealing 
time from 1 to 15 h always resulted in a significant sharpening of the DTA 



Chromel*alumel differential thermocouple voltage (pV) 
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curves in the vicinity of the apparent solidus temperatures, which indicates 
that the annealing effected alloy homogenization by diffusion. On the other 
hand, the changes caused by Increasing the annealing time beyond 30 h were 
usually negligible, as shown for samples with x • 0.1 and 0.2 in Figures 4 
and 5, respectively, where and correspond respectively to the annealing 
and liquidus temperatures, and the symbols 1, 2, 3, and 4 designate the DTA 
curves that weie obtained following annealing times of 1, 15, 30 and 60 h, 
respectively. 

DTA curves for fully annealed samples with x » 0.1, 0.2, 0.3, 0.4, 0.6, 
0.7, 0.8, and 0.9 are shown in Figure 6, and the collection of melting curves 
is shown in Figure 7. The horizontal axis shows the measured Sb-reference 
ampule temperature, T, and the vertical axis gives differential temperatures, 
S, for the alloy ampule as measured with respect to T. 



Figure 4. Differential thermal analytit data for Hgg gCdQ ‘]Te. 
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Figure 6. DTA data for x » 0.1, 0.2, 0.3, 0.4, 0.6, 0.7, 0.8 and 0.9 Hg.j,^Cd^Te niloyi. 
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Figure 7. Hg|_j^CdjjTe melting curves as observed. 

Considering radiation heat transfer only, the heat balance at the outer 
wall of the sample ampule yields 


e T " - e T,^ (e T " - e O. 

FF pi r^pl ss'^ 

^ q 


( 1 ) 


where Tp, Tp and Tg are the furnace, sample thermocouple, and sample 

temperatures, respectively, e , e , and e are the furnace, sample, and quart: 

F s p 

emissivltles, respectively, and r_ and r_ are the inner and outer radii, 

5 q 

respectively, of the quartz ampules. The measured differential thermocouple 
signal is given by S ■ Tj^ - T. For steady-state conditions where the furnace 
temperature increases at a constant rate, T ■ Tp + P' , where P' is a 
constant. Using the approximations 


Tj = (t + S]^ - + 4ST^ 


( 2 ) 


and 


= (t + P']^ - + 4P'T^, 

r 


(3) 
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Equation (1) yields 
™ A 



r/ r \ 

e 


« T^ 

(1 +-^1 

_E 

_ -a 



e 

r e 


i\ 9/ 

s 

s s 


+ AST' 




e , r e_ 

_E . 4p-t 3 _a _F _ 
e re 

s s s 


(A) 


A first-order expansion In f/T and P'/T yields 




(5) 


where a,, ao, and a-i are functions of r /r_, e /e^, and c„/c . 

j qs ps ra 

For HgTe the observed empirical correlation between T, S, and T_ is given by 

8 


T = T + 2S - P. (6) 

s 

For Cdle also, S2 was found to be 2. In the case of a nonzero baseline, S^, 

T„ Is given by 
8 

T - T - S + 2S - P. (7) 

s o 

P is primarily the temperature drop across the ampule wall, but this term also 
Includes thermocouple calibration corrections. In the absence of a phase 
transition, ■ T + - P, l.e., the sample temperature is equal to the 

reference thermocouple temperature plus the differential temperature minus the 
temperature drop across the ampule wall. 

For each sample, P was determined by the difference between the observed 
and actual melting points of the Sb-reference. The values measured for P 
using other calibration samples with known melting points ranged from 1.3°C to 
3.A°C, which agree reasonably with finite-difference calculations that include 
both radiation and conduction heat transfer. 

The solidus temperatures were determined by replottlng the data of 
Figure 6 as a function of Tg Instead of T, using Equation (7). The solidus 
temperatures for 0.1 < x < 0.7 were determined by extrapolating the initial, 
straight-line part of the replotted melting curves to the extended baseline. 
For X - 0.8 and 0.9, the initial departure from the baseline was used because 
the melting curves had no straight-line sections even when replotted as 
functions of T_. The liquldus points, in all cases, were considered to be the 
temperatures at which the melting curves broke suddenly toward the baseline. 
The observed phase equlllbrlun data are listed in Table 3. 
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TABLE 3. MEASURED PHASE EQUILIBRIA DATA FOR Hg^.^Cd^Te. 


Composition 

(x) 

So 

at solidus 
(«C) 

So 

at liquidus 
l°C) 

P 

l°C) 

Observed 

solidus 

temperature 

(®C) 

Observed 

liquidus 

temperature 

("O 

0.1 

0.5 

-1.5 

2.9 

690 

733 

0.2 

0.4 

-1.0 

2.9 

706 

791 

0.3 

1.4 

0.7 

2.9 

727 

841 

0.4 

0.9 

0.3 

3.4 

748 

882 

0.6 

0.3 

-0.3 

1.3 

810 

952 

0.7 

0.7 

-0.4 

2.9 

850 

990 

0.8 

6.9 

4.6 

2.9 

904 

1022 

0.9 

3.4 

-1.0 

2.9 

950 

1051 


aP03-07S0'4 


Each ampul- contained a small free volume. Because the partial vapor 
pressure of Hg over the alloy melt is much larger than the partial pressures 
of Cd and the preferential evaporation of mercury into the ampule free 

volume alters both the Te/metal and Cd/Hg fractions. As suggested by the 
decrease of the binary Hg-Te and Cd-Te liquidus temperatures with increasing 
Te composition immediately above the 50 at.% Te composition, an 
increase in the Te/metal fraction can significantly affect the liquidus 
temperature. The Cd-Te liquidus temperature decreases at *»10°C/at.% Te^»^» 
and the Hg-Te liquidus temperature decreases at 1.2°C/at.% Te.^’^^ 

The following procedure was used to correct for the small Te-metal excess 
in the melt caused by preferential evaporation of Hg into the ampule free 
volume. Following the casting of an alloy, the free volume was estimated by 
measuring the length of the free space in the ampule. The amount of Hg in the 
gas phase was determined by using Stelninger's^ empirical relationship for the 
Hg partial pressure, P, given by 

P = exp(21.732 - T^/T), (8) 

where P^ = 1 Pa and T^ = 7.149 K. The number of moles of Hg, nj^g, remaining 
in the melt was obtained by subtracting the amount of Hg in the vapor from the 
amount loaded into the ampule. The actual Te/metal mole fraction in the melt 
at the observed liquidus temperature is then n^g/(njjg + n^^j + n.j.g) , where n^.^ 
and njg are respectively the number of moles of Cd and Te in the melt. The 
relatively small amounts of Cd and Te 2 in the vapor were neglected, and the 
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values of and Oj-g were assumed to be the same as originally placed in the 
ampule. The corrections to ti;e observed llquldus temperatures were then 
obtained from the expression 


“Te 

"llg "Cd "Xe_ 


-0.5 (1.2 + 8.8x) , 


(9) 


which is based on a linear interpolation between Hg-Te and Cd-Te of the rate 
of change of the llquldus temperature with respect to the Te concentration for 
Te mole fractions slightly above 0.5. The measured free volumes, the 
calculated atomic Te concentrations, the estimated temperature corrections for 
free-volume effects, and the corrected llquldus temperatures are summarized in 
Table 4. 

Table 5 lists the llquldus and solidus temperatures determined in this 
Investigation. Table 5 also gives the corrected x-values based on 
X ■ *'cd'^^"Cd "Hg^* phase diagram data are plotted in Figure 8. 


TABLE 4. FREE-VOLUME CORRECTION PARAMETERS AND CORRECTED 
LIQUIDUS TEMPERATURES. 


composition 

(x) 

Capsule 

free 

volume 

(cm^) 

Te fraction 
(%) 

Liquidus 

temperature 

correction 

(°C) 

Observed 

liquidus 

temperature 

i°C) 

Corrected 

liquidus 

temperature 

m 

0.1 

0.19 

50.02 

0 

733 

733 

0.2 

0.26 

50.04 

0.1 

791 

792 

0.3 

0.69 

50.11 

0.4 

841 

841 

0.4 

0.88 

50.21 

1.0 

882 

884 

0.6 

0.96 

50.33 

2.1 

952 

954 

0.7 

0.74 

50.27 

2.0 

990 

992 

0.8 

1.46 

50.31 

2.5 

1022 

1025 

0.9 

0.87 

50.20 

1.9 

1051 

1063 


0NMn»l 


16 



( 4 ) 


ORIGINAL PAGt 
OF POOR QUALITY 

TABLE 5. LIQUIDUS AND SOLI DUS TEMPERATURES FOR Hg^.^CcIxTe. 


Nominal 

composition 

(x) 

Actual composition 
at solidus 
(x) 

Solidus 

temperature 

(°C) 

Actual composition 
at liquidus 
(x) 

Liquidus 

temperature 

l“C) 

0.1 

0.100 

690 

0.100 

733 

0.2 

0.200 

706 

0.200 

792 

0.3 

0.301 

727 

0.302 

841 

0.4 

0.402 

748 

0.404 

Cn/i 

C.6 

0.604 

810 

0.607 

954 

0.7 

0.704 

850 

0.707 

99' 

0.8 

0.806 

904 

0.810 

io; 

0.9 

0.905 

950 

0.907 
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Figures. Hg^.^Cd^Te phaie diagram. 
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2.4 Error Analysis and Comparison vdth Previous Results 


The major recognized sources of error in the phase diagram data are the 
thermocouple calibrations, the subjectivity In selecting the proper phase- 
transformation points from the DTA curves, and the uncertainties In the ampule 
free-volume corrections. 

Because the Sb melting point (630. 7°C) was user", to calibrate the 
thermocouples during each DTA run, the calibration should be quite accurate 
near that temperature and decrease In accuracy with Increasing temperature. A 
type-K thermocouple calibrated at 630. 7°C should be accurate to within ±2°C up 
to 930°C.^^‘^^ The measured melting temperature of a Ag standard sample 
differed from the IPTS-68 value of 930°C by not more than ±3®C in tests of 
three different .hermocouple combinations used during the courra of this 
Investlgal on. Therefore, the expression 3(T - 63i°C) / (930-63n was used as 
the calibration uncertainty at temperature T, which Is equivalent to a\ 
uncertainty at the Ag melting point of 3.3°C and Is considered adequate to 
Include any graph-paper calibration errors as well as thermocouple calibration 
errors. For each sample, the calculated uncertainties arising from the 
temperature measurement uncertainties are listed In Table 6. 


TABLE 6. UNCERTAINTIES OF THE SOLIDUS AND LIQUIDUS TEMPERATURES. 


Composition 

(x) 

Solidus temperatures 


Liquidus temperatures 


Temperature 

calibration 

(®C) 

Reading 

DTA 

curve 

(OC) 

Total 

(°C) 

Temperature 

calibration 

(°C) 

Reading 

DTA 

curve 

<0C) 

Free 

volume 

effect 

(“Cl 

Total 

(®C) 

0.1 

0.6 

1.0 

1.6 

1.0 

3.0 

0 

4.0 

0.2 

0.8 

1.0 

1.8 

1.6 

1.0 

0 

2.6 

0.3 

1.0 

1.0 

2.0 

2.1 

1.0 

0.1 

3.2 

0.4 

1.2 

1.0 

2.2 

2.5 

1.0 

0.3 

3.8 

0.6 

1.8 

2.0 

3.8 

3.2 

hO 

0.5 

4.8 

0.7 

2.2 

1.0 

3.2 

3.6 

1.0 

0.5 

5.1 

0.8 

27 

2.0 

4.7 

3.9 

2.0 

06 

6.6 

0.9 

3.2 

2.0 

1 / 

4.2 

2.0 

0.5 

67 
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The second source of error is in determination of the critical points of 
the DIA curves. For the solidus points, the extrapolation of the S vs Tg 
curve to the baseline for alloys with 0.1 < x < 0.7 and the choice of the 
onset of deviation from the baseline for alloys with x = 0.8 or 0.9 are 
somewhat subjective. Similarly, there were uncertainties in identifying the 
liquidus points. The latter uncertainties were about 1°C except for the 
X = 0.1 sample, for which the liquidus onset occurred where the curve had a 
very steep slope, and for the x = 0.8 and 0.9 samples, for which the sharpness 
of the liquidus onset was somewhat more rounded as shown in Figure 6. 

Estimates of subjectivity in reading the DTA curves are given in Table 6 under 
the columns headed "Reading DTA Curve." 

The uncertainties listed for the free— volume corrections are 25% of the 
total free- olume corrections. 

T>ie total uncertainties listed in Table 6 are the sums of the individual 
uncertainties and are believed to represent the maximum possible errors in the 
reported phase-equilibria temperatures. 

The solidus and liquidus temperatures are compared with the results of 
previous investigations in Figure 9. The data of Ray and Spencer^ ^ and Blair 
and Newnham^® were derived from DTA measurements. Harman's^' data are from 
chemical analyses of first-to-freeze regions of molten alloy samples. 
Steinlnger's^ data were derived from thermal analyses of samples in a high- 
pressure apparatus. 

A DTA measurement tends to yield too low a solidus temperature for a 
given X because of possible Inhoraogeneities in the sample. The inherent error 
of the first-to-freeze method tends to yield too low an x value for a given 
temperature. Therefore, on a temperature-composition plot, a solidus curve 
determined by the first-to-freeze compo itions should fall at or above the 
true curve, and a solidus curve determined by DTA measurements should fall at 
or below the true curve. The present DTA data establish a solidus curve that 
nearly coincides with the lower limit of Harman's first-to-freeze data, and 
the two independent sets of data thus are corroborative. 
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Figure 9. Comparison of present phase-equilibria data with previous results. 


For most compositions, the solidus and liquidus temperatures determined 
by this investigation are much larger than those reported by Ray and 
Spencer. A possible explanation for the descrepancies is advanced in 
Reference 20. 


2.5 Liquid-Solid Equilibrium Parameters 

The equilibria temperatures in Table 5 were used to calculate the 
temperature and composition dependence of the liquid/solid interface 
distribution coefficient, k(T) = x^(T)/x^(T), where x^ and x^ are the solidus 
and liquidus compositions, respectively, fo facilitate the calculations, 
analytical expressions were developed for x^(T)/x^(T). Reasonably good fits 
to the phase diagram data were given by the following functional forms: 


x^(T) 


C. slnl-r 


/tt „ . A „^l/2\. „ , /9T - 5618°c\ 


C.T* 


1/2 


(IJ) 


and 
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x^(T) = DjT* + 02 !*^ + D^T*^ + 

where T* = (T - 67C°C)/412°C and the and are fitting constants. The 
values of the constants giving the best fits are given in Table 7. The solid 
curves drawn in Figures 8 ana 9 are plots of Equations (10) and (11)* The 
calculated values of k(T) as functions of temperature and x^(T) are shown, 
respectively, in Figures 10 and 11. The failure of k(T) to go through the 
point (T, k) = (670°C, 1) or the point (x^, k) = (0, 1) is an artifact of the 
solidus fit shown in Figure 8. 

TABLE 7. BEST VALUES FOR THE FITTING CONSTANTS Cj AND Dj. 


i Cj Dj 

1 0.380857 0.607640 

2 0.086277 0.077209 

3 1.200139 0.696’ 67 

4 -0.665231 -0.381683 


GP0M7S04 



Figure 10. Temperature dependence of the interface distribution coefficient 
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Figure 1 1. Composition dependence of the distribution coefficient 


To compare the observed DTA signals vd.th those expected from the alloy 
properties, the equilibrium liquid fraction, Ml(T), was determined by using 
the lever rule, Ml(T) * [x^(T) - Xq]/[x^(T) - x^(T)], and the analytical 
expressions of Equations (10) and (11) for x^ and x^. The temperature and 
composition dependences of Ml(T) are shown in Figure 12, and their temperature 
derivatives are shown in Figure 13. Figure 14 shows the observed normalized 
DTA signals for the various alloy compositions. A comparison of the results 
shown in Figure 14 with the results in Figure 13 shows a one-to-one 
correlation between the relative magnitudes of th DTA signals and "he 
variation of dM^^/dT with composition. In principle, it the heat tran.'^fer 
characteristics of the furnace/sanple system are well knowr», the compositional 
and temperature variations of dM^/dT and thus Ml(x^. x^) can be calculated 
from the measured thermal arrest curves and vice versa. Approximate heat 
transfer calculations were performed, and in general, Che calculated and 
observed DTA curves for the various alloy compositions were comparable. 
However, the detailed features of the calculated curves were sensitive to 
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small variations in the values of the heat transfer parameters, some of which 
are not precisely known. Meaningful decon»folutions of the measured DTA curves 
would require a precise thermal model for the sample/furnace system. 



Temperature (°C) 


OFO3-O7S0-22 

Figure 12. Temperature dependence of the alloy liquid fraction for various alloy compositions. 
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3. PHASE EQUILIBRIA CALCULATIONS 


Iterative phase-equilibria calculations, based on the regular associated 
solution (R.A.S.) theory^*^ were performed, and a set of R.A.S. parameters was 
obtained by simultaneously least-squares-f itting the binary^ and 
pseudobinary phase diagrams. Using the set of R.A.S. parameters, the 
activities and partial pressures of Hg, Cd, and Te vapor over pseudobinary 
melts were calculated. 


3.1 Theoretical Model 

The term, regular solution, was used by Hildebrand^ ^ to describe mixtures 
whose behavior was found experimentally to have certain regularities, and the 
definition of the term is essentially empirical. Subsequently, Guggenheim 
used the term to describe any mixture of molecules for which the entropy is 
that of an ideal solution but for which the interchange energy is not zero. 

The interchange energy is related to the nearest-neighbor pairwise interaction 
energies between atoms. Subsequently, Jordan^ considered the consequences of 
the association of type A and type B atoms In the liquid, i.e., the existence 
of AB species, and developed a thermodynamic description for such liquid 
solutions. He termed his formulation the regular associated solution (R.A.S.) 
theory. We assumed that the binary alloys Hgj_j^Tex hence 

the ternary Hg-Cd-Te alloys, are regular associated solutions, and used 
Jordan's formalism, as extended to ternary alloy systems by Szapiro,^ to 
determine the pertinent thermodynamic phase-equilibria parameters for the 
alloy system. 

In the context of the R.A.S. theory, the liquidus and solidus ternary 
equations for an ABC ternary system are given by^: 


[1 - y]exp 


W 2 
RT ^ 


AC 

X * 

AC 


|-AHi 
X exp 


AH, + 00 X 2(1 “ 2x^) + ~ 0*5oj^) (1 - 2x^) AS 


RT 


AS^- 
" “r“ ’ 


( 12 ) 
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y exp 


[I? 


X ^ 

BC 


X exp 


■AH 2 + tox^ (1 - 2 X 2 ) + + 012*2 ~ 0.5a2> (1 - 2x^) AS 


RT 


R J ’ 


( 13 ) 


where: 

y 

xj, X2» X3 

*AC» *BC 
T 

AS ^ y A S 2 
AHj, AH2 
R 

“1 

“2 

u 

and 

W 


mole fraction of complex BC In the solid phase, 
atom fractions of A, B, and C, respectively, 
mole fractions of AC and BC in the liquid phase, 
absolute temperature, 

entropies of fusion of the binary systems AC and BC, 
heats of fusion of the binary systems AC and BC, 
the universal gas constant, 

the interaction parameter for the liquid binary 
system AC, 

the interaction parameter for the liquid binary 
system BC, 

the interaction parameter between species A and B, 

the interaction parameter between species AC and BC 
in the solid. 


Asterisks denote the standard states of the liquid phases (stoichiometric 
compositions in the binary systems AC and BC). The Interaction 
parameters a^, u, and W are related to the interchange energies for the 
various pairwise interactions. 


For the binary case, y “ 0, X 2 =• 0, and xj + X 3 ■= 1. Letting x = X 3 , 

Equation (12) for this case reduces to 


exp 


AHj + 2o (x - 0.5)^ 
_ 



(14) 
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where 


1 - [1 - 4(1 - x)x/q + 

1/2 * 

1 + [1 - 4(1 - x)x/(l + tc)l ' 


( 15 ) 


X 


* 

AC 


AC 


(x “ 0.5) , 


(16) 


and If Is given in terms of the temperature- independent equilibrium constant, 
K, by 


/ s' \ 

N iv^r 


The parameter B is defined by Equation (17) and is a measure of the degree of 
association in the liquid. For a completely associated solution, ic * 0, and 
thus B « 0. When the solution is completely dissociated, ic - » and B • The 
activities of the species a^^ and a^ for the binary melt are given by 



(18) 


(19) 


where y, and y„ are the activity coefficients, 
A C 


[1 + 


”ac' - 


‘AC' 


( 20 ) 


and 


^C ” ^3 *AC^ ■ ’'AC* 


( 21 ) 


For the pseudobinary cut, x^ ■ 0.5 = x^ + x^, so that Equations (12) and (13) 
become 


(1 - y) exp 


W 2 * 

’'ac 

"AHj + (0X2(1 “ 2 xj^) 

ASj" 

RT ^ 

” X * 

’'ac 

FT 



( 22 ) 
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and 


[PT 


y exp — (1 - y)' 


^BC 


exp 


rAH2 +<DXj(l - 2X2) ^^2! 

[ RT FJ ’ 


(23) 


where 


AC 


Xj^ x^ (1 + o) 

"1 + 


(24' 


*BC 


X2 Xj, (1 + 0) 


(25) 


and a (= x^^ + x^^) is the solution of the cubic equation 


+ - 


+ 0.5 (1 - 0 ) ’ <2+0.5 (1 - a) _ ^2^ ‘ 


(26) 


The parameters <j(AC) and < 2 (BC) are defined by Equation (17), and the 

parameters x* (i = A, B) are defined by Equation (16). 
l.C> 

The activities a^^, ag, and a^ for the pseudobinary melt are given by the 
expressions: 


*A ■ Va ■ + ’‘AC ’‘bc> ' ■‘acJ’a' 


= XV, 


[x-(l + X. . + x_) - x_.]y-. 


and 


‘C ■ »C ■'c ■ - -AC - 


where y^, and y^ are the activity coefficients given by 


( 0 X 2 + 0.25 + (o) + - O 2 ) (0.5 X2> 



(27) 

(28) 


(29) 


(30) 


( 31 ) 




and 
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exp 


2 2 

“2^2 ®2 ~ ^1*2 

RT 


(32) 


3.2 Regular Associated Solution Parameters for the HgCdTe Alloy System 

The R.A.S. parameters a^, 02 * < 2 » <*>» and V? for the HgCdTe alloy 

system were determined by simultaneously least-squares-f ittlng the metal-rlch 
segments of the binary liquldus curves^”^ and the pseudobinary (Figure 8) 
phase diagram. To facilitate the computations. Equations (22) and (23) were 
rewritten in the following forms; 


+ 00X2 (1 - 2 x^) - Wy^ 


ln(l 


. 1 AC 

y) - In ^ 

"'ac 


+ AS 


1 


and 


(33) 


AH 2 + o)x^ (1 - 2 X 2 ) - W (1 - y)2 


R 

In y - In 

+ AS- 



2 


(34) 


These equations and Equation (14) were then fitted to the experimental data to 
obtain the best values of the parameters by minimizing the sum: 


S 


CdTe 


(v‘> 


CdTe\^ 
-T^(e)j + 


' HgTe HgTe\^ 
T^(c) - T^(e) 1 


+ 


(t 


HgCdTe 

31 


HgCdTeX 
T,(e) j 


HgCdTe 

21 


HgCdTeX 

Ti(e) j 


(35) 
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The symbols c and e designate, respectively, the calculated and experimental 
values* 

Table 8 lists the calculated values of the parameters as well as the 
pertinent thermodynamic quantities used in the calculations* 

Table 9 indicates the relative accuracy of the fit to the measured 
pseudobinary data* The calculated mole fractions of the various species and 
complexes as functions of the solidus temperature are given in Table 10 to 
Illustrate the relative effects of association on the Hgj_jjCdjjTe solutions* 


TABLE 8. R PARAMETERS AND THERMODYNAMIC DATA FOR THE Hg-Cd-Ta 
SYSTEM. 



R.A.S. parameten 

Thermodynamic data 

-0.9716 

0 (^ - 953.0 cal/mol - 3.987 IcJ/mol 
K2- 0.1201 

02 ■ 9051 cal/mol - 37.87 kJ/mol 
(u - 971 .7 cal/mol - 4.066 kJ/moi 

W - 836.0 cal/mol - 3.498 kJ/mol 

ASy - 9.2 e.u. - 38.5 J/mol 
• 8.9 e.u. - 37.2 J/mol 
- 8675.6 cal/mol = 36.299 k J/mol 
AH 2 - 1 201 2 cal/mol - 50.258 kJ/mol 
Tf^ -943K<a> 

Tf -1365K<®> 


**^A. Laugiar, Ravu* Phytiqu* AppliquM 8, 2S9 (1973). QP0$47S0-> 

tubicrjpt 1 rtftrs to HgTe and tha lubtcript 2 refers to CdTe. 
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TABLE 9. A COMPARISON OF EXPERIMENTAL 
AND CALCULATED TEMPERATURES 
FOR HQi.^CdxTe. 


^exp 

IK) 

^2 

(K) 

T3^ 

(K) 

973 

968 

981 

998 

991 

1003 

1023 

1012 

1026 

1048 

1032 

1052 

1073 

1053 

1074 

1098 

1074 

1098 

1123 

1097 

1122 

1148 

1121 

1145 

1173 

1149 

1168 

1198 

1180 

1192 

1223 

1213 

1217 

1248 

1254 

1242 

1273 

1298 

1268 

1298 

1358 

1294 


The quantttiet deiignated by t were calculated 
using Equations (33) and 1341 ; the parameters are given 
in Table 8. 
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TABLE 10. MOLE FRACTIONS OF SPECIES IN LIQUID 
PHASE AS FUNCTIONS OF LIQUIDUS 
TEMPERATURE. 


(K, *Hg *Cd *Te *HgTe xCdT* 


973 

998 

1023 

1048 

1073 

1098 

1123 

1148 

1173 

1198 

1223 

1248 

1273 

1298 


0.387 

0.379 

0.370 

0.359 

0.348 
0 335 
0.319 
0.302 

0.282 

0.258 

0.231 

0.200 

0.164 

0.121 


0.001 

0.002 

0.003 

0.004 

0.005 

0.007 

0.010 

0.013 

0.017 

0.022 

0.029 

0.038 

0.050 

0.066 


0.388 

0.380 

0.373 

0.363 

0.354 

0.3^2 

0.329 

0.315 

0.299 

0.281 

0.260 

0.238 

0.214 

0.187 


0.198 

0.189 

0.179 

0.169 

0.158 

0.147 

0.134 

0.120 

0.106 

0.091 

0.075 

0.059 

0.044 

0.028 


0.026 

0.050 

0.075 

0.104 

0.134 

0.169 

0.207 

0.250 

0.296 

0.348 

0.405 

0.465 

0.529 

0.597 
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3.3 Component Partial Pressures over the Pseudoblnarv Melt 


Following Jordan, the partial pressures of the components over the 
pseudobinary melt can be related to the activities by the expressions 



Hg 



( 36 ) 


Cd 


®Cd ^Cd’ 


(37) 


and 


P 


Te 


2 




(38) 


where P°^, P°^, and P°^ are the vapor pressures over the standard states and 
the respective activities a^g, acd» ®Te given by Equations 
(27-29). Using Equations (36-38) and the parameters of Table 8, the component 
partial pressures were calcula*^ ?d as functions of the llquldus temperature, 
and the results are summarized In Table 11 and Figure 15. Figure 16 
Illustrates the temperature variation of the Hg partial pressure over the 
pseudobinary melt as a function of the Cd fraction. The calculated Hg partial 
pressures are In reasonable agreement with Stelnlnger's experimental 
results.^ However, the validity of the calculated R.A.S. p3r->’;,i.t .i values and 
thus the calculated partial pressures are sensitively dependent upon 
assumptions concerning the melt composition. For example, the derivation of 
Equations (12) and (13) assumes negligible Hg-Hg, Cd-Cd, and Te-Te 
associations In the melt. The validity of at least the last of the 
assumptions Is highly suspect. 
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TABLE 11. CALCULATED COMPONENT PARTIAL PRESSURES OVER THE 
H9^^Cdj^T« MELT AT THE LIQUIDUS TEMPERATURE. 



Composition ^Hg 
(K) (x) 

«Cd 

»Te 

po (a) 
(atm) 

pO (b) 

(atm) 

(atm) 

’’Hg 

(atm) 

•*Cd 

(atm) 

'*Te2 

(atm) 


973.2 

0 044 

0.421 

0.009 

0.457 

50.6 

0.457 

0.039 

21 . 3 ( 17 . 5 ) 

0.004 

0.008 

998.2 

0.084 

0 396 

0.016 

0.463 

60.8 

0.626 

0.056 

24.1 ( 21 . 0 ) 

0.010 

0.012 

1023.2 

0.124 

0.374 

0.023 

0.468 

72.5 

0.844 

0.078 

27.1 ( 25 . 0 ) 

0.020 

0.017 

1048.2 

0.170 

0.350 

0.031 

0.473 

85.6 

1.120 

0.108 

29.9 (29 

0.035 

0.024 

1073.2 

0 . 21G 

0.327 

0.039 

0.476 

100.0 

1.470 

0.146 

32.8 ( 34 . 5 ) 

0.057 

0.033 

1098.2 

0.268 

0.302 

0.048 

0.479 

117,0 

1.910 

0.196 

35.3 ( 40 . 4 ) 

0.091 

o.cir 

1123.2 

0.324 

0.278 

0.057 

0.480 

135.0 

2.440 

0.259 

37.5 ( 47 . 0 ) 

0.140 

0 . 06t 

1148.2 

0.384 

0.252 

0.068 

0.479 

155.0 

3.090 

0.338 

39.1 ( 53 . 5 ) 

0.210 

0.078 

1173.2 

0.446 

0.227 

0.080 

0.475 

177.0 

3.880 

0.436 

40.1 ( 61 . 0 ) 

0.309 

0.098 

1198.2 

0.514 

0.200 

0.094 

0.468 

201.0 

4.820 

0.556 

40.2 ( 69 . 0 ) 

0.451 

0.122 

1223.2 

0.586 

0.172 

0.110 

0.456 

227.0 

5.930 

0.703 

39.0 (' 0 ) 

0.654 

0.146 

1248.2 

0 660 

0.143 

0 . 1 3C 

0.439 

255.0 

7.250 

0.880 

36.5 

0 943 

0.170 

1273.2 

0.736 

0.113 

0.155 

0.415 

286.0 

8.780 

1.090 

32.3 

1.360 

0.188 

1298.2 

0.816 

0.081 

0.188 

0.385 

318.0 

10.600 

1.340 

25.7 

1.980 

0,199 


Values in parentheses correspond to Steimnger’s experimental results ^Ref 7) 
^Hg* ^Te^ vapor pressures over the elemental melts 


(a) Pu was obtained from the equation. loq,rtP?i <atm) = 4.89'i — T7TT7. which gives with an accuracy of better than 

rig ' I U Mg T I K ) rig 

1% [S. Sugawara et al.. Bull, Jap. Soc. Mech. Eng. 5, 711 (1962)1 . 

(b) P^ was obtained from the equation. (atm) - 5.106 - which gives P^^ with an accuracy of better than 

0.5% (D.R. St’jll and G.E. Sinke, Thermodynamic Properties of the Elements , (Am. Chem. Soc., Washington, O.C. 1956)1 

(c) Pje- was obtained from the equation, log-^P^ (atm) ^ 4 7191 - ^ (An. N. Nesmeyanov, Vapor Pressure of 

A 1 U iCj T 1 1C 1 

the Elements (Academic Press. New York. 196311. GP030750-12 
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Figure 15. Component partial pressures over the Hg^.^Cd^Te melt. 
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Figure 16. The temperature variation of Hg partial pressure over the Hg^.j^Cd^Te melt 
at a function of Cd fraction. 
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4. HIGH-TEMPERATURE-CRADIENT DIRECTIONAL SOLIDIFICATION 


4.1 Alloy Preparation 

alloys were prepared by reacting 99.9999% pure elemental 
constituents In sealed, evacuated, fused-quartz tubing. The tubing dimensions 
were 5-mm l.d. x 10-mm o.d. To prepare an ampule, the end of the tubing where 
solidification was to begin was tapered to a point to enhance the probability 
of single-crystal growth. For the initial experiments, a 5-nmi l.d. x 10-mm 
o.d. quartz tube was fused to the tapered end of the ampule to serve as a 
coupling sleeve to a 5-nnn diam quartz-rod pedestal, which supported the ampule 
during the growth process. This coupling arrangement was susceptible to 
mechanical failure, and after several failures, a configuration was used in 
which a 4-nm diameter silica rod is fused directly to the tapered end of the 
ampule. The solid silica rod inserts into a 5-nmi l.d. silica tube, which 
serves as a pedestal. 

Following fabrication, the ampules were cleaned thoroughly, with the 
final cleaning being done with a dilute HF solution, evacuated, and baked at 
1150°C. The ampules were loaded with precisely weighed amounts of Hg, Cd, and 
Te. Following loading, the ampules were sealed and the alloys were cast as 
described in Section 2.1. 

Altogether, 12 crystal-growth ampules (L1-L12) containing alloys with x = 
0.2 were prepared. The ampules L1-L5 contained stoichiometric proportions of 
the constituent elements. The alloy preparations yielded highly p-type ingots 
because of the Te excess in the melt caused by the preferential Hg evaporation 
into the ampule free volume. For ampules L6-L12, excess mercury was added to 
compensate for the preferential mercury evaporation. Ampules LI, L2, and L8 
failed during casting. Ampules L9, LIO, and Lll showed wetting of the ampule 
wall by the alloy, and they were discarded. 

4.2 Experimental Arrangement and Procedures for Unidirectional Crystal Growth 

The precast alloys were regrown by unidirectional solidification In the 
growth apparatus illustrated in Figure 17. The ampule is mounted, with its 
tapered end down, atop a quartz pedestal and remains stationary during the 
growth process. Two reslstlvely heated, tubular furnaces provide for the 
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Figure 17. Bridgman crystal-growth furnace assembly. 

temperature control of the hot and cold zones* The furnace assembly travels 
upward on vertical ways so that the crystal grows progressively upward from 
the bottom, tapered end of the ampule. Each furnace Is provided with a heat- 
pipe isothermal furnace liner to provide two well-defined heated zones. The 
hot-zone heat-pipe is sodium-charged, and the cold-zone heat-pipe is either 
sodium or mercury charged depending on the temperature requirements. The 
heat-pipes are separated by a heat barrier. Besides providing the necessary 
thermal isolation between the heated zones, the barrier also provides a semi- 
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adiabatic region. The barrier consists of 4-mm thick zirconia-f elt or various 
thicknesses of fire-brick disks sandwiched between two 25-wm thick polished 
platinum foils on the hot side and two 25-um aluminum foils on the cold 
side. The diameter of the circular opening In the heat barrier Is 12-mm to 
accept 10-mm diameter ampules. 

Prior to a given growth run, a dummy ampule with an attached thermocouple 
Is placed in the upper zone of the furnace, and the two heated-zone 
temperatures are raised to their desired values. The furnace assembly Is 
translated upward, and the thermocouple output is recorded as a function of 
furnace position to establish the temperature variation along the furnace 
assembly axis. The furnaces are then cooled to ambient tec^erature, the 
growth aiq)ule Is positioned In the upper zone, and the zone temperatures are 
increased to the previous values. Prior to growth Initiation, the molten 
alloy Is maintained at temperature for 24 h to assure melt homogenization. 

4.3 Bridgman Crystal Grotrth Runs 

Crystal growth runs were completed for ampules L3, L4, L5, L6, and L7 at 
furnace translation rates of 0.597, 1.12, 5.62, 0.310, and 0.068 vm/s, 
respectively. The longitudinal temperature profiles in the furnace at the 
initiation of growth of Ingots L3 and L4 are shown In Figures 18 and 19, 
respectively, and the furnace profile for ingots L6 and L7 is shown in Figure 
20 . 

The recrystallized Ingots were removed from the ampules by dissolving the 
quartz in concentrated HF. Sample L5 showed excessive wetting between the 
ampule and the Ingot and was not evaluated further. The remaining ingots were 
sliced into wafers by a wire saw for characterization. 
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Figure 18. Temperature profile of crystal-growth furnace during 
growth of ingot L3. 
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apo3-arM'4i 

Figure 19. Temperature profile of cryatal-growth furnace during growth of ingot L4. 
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Figure 20. Temperature profile of crystal -growth furnace during 
growth of ingots L6 and L7. 


40 



4.4 Ingot Characterization 


4.4.1 Determination of Alloy Composition by Mass Density Measurements 

The crystal lattice spacing of the Hgj_j^CdjjTe alloys varies so little 

with X that its measurement cannot be used readily to calculate precise values 

for the compositions of specimens. For this study, the average x-values of 

slices transverse to the growth axis were calculated from the measured mass 

densities and values of the crystal lattice constant published by Woolley and 
'>4 

Ray. 

The lattice constant data of Woolley and Ray can be represented by the 
cubic equation 


a^ = 0.646153 + (0.73671 x + 1.90501 x^ - 0.69347 x^)/1000. 


(39) 


where a^ is the cubic lattice constant in nanometers and x is the mole 
fraction of CdTe in the Hg 2 _ 3 ^CdjjTe alloy. There are, on the average, 4x Cd 
acorns, 4(l-x) Hg atoms, and 4 Te atoms in each cubic unit cell. Thus, the 
lattice constant is related to the mass density, p^, by the expression 




1/3 


Hg 


Te' 


(40) 


where Mjjg, and are the atomic masses of Cd, Hg, and Te, 

respectively. For a given density, the value of x can be determined by 
simultaneously solving Equations (39) and (40). 

The mass density of a specimen was determined by obtaining its weight, 
Wg, in air and its weight, when it is immersed in water. The specimen 
density is given by the relation 


W p - W p 
- a w w a 

^m " W - W ’ 


(41) 


where p is the mass density of water and p is the mass density jf air. 

W 8 
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Near x * 0.2, an uncertainty of 0.0022 in leads to an 

uncertainty of about 0.001 in calculating x. The uncertainty in can be 
approximated by 


Ap 


m 


* ‘’m “K * “g 

w - w 

a w 


( 42 ) 


where AW , AW , and Ap are the uncertainties in W , W , and p , 
a w w a w w 

respectively. Because AW and Ap are usually much smaller than AW , 

a w 


Ap. 


m 


p AW 
m w 

W - W 
a w 


2 

p AW 
m w 

W p 
a w 


p^ AW 
m w 

W 


(43) 


The uncertainty in x. Ax, is thus given by 


Ax 


AW p^ 
w m 

2.2 W 


» 


(44) 


where the density is in the units of g/cm^. For specimens that weigh about 
0.5 g, AW^ must be less than 0.00002 g to give x to within ±0.001. 

The method described by Bowman and Schooner^^ was used for the mass 
density determinations. The Hg^.j^Cd^Te specimen is immersed in water, and the 
water is boiled briefly to eliminate gases dissolved in the water and to rid 
the specimen of minute air bubbles that may be sorbed. After the sample and 
the water have cooled to room temperature, the weight of the immersed specimen 
is measured with a balance that weighs with a precision of ±0.000001 g. 
However, the precision of the weight determination in water is much poorer 
than this because of the surface-tension forces upon the 0.013-mm dlam wire 
that is used to suspend the specimen from the balance arm. Some of the 
suspension wires that were prepared permitted the weight in water to be 
determined with an accuracy of ±0.00001 g, but usually an accuracy of only 
±0.00002 g is achieved. To do significantly better requires an Inordinate 
amount of time. 
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The following expression, from Reference 2L-, was used to calculate the 
density of the water In which the sample Is Immersed: 


1 - 


(t - 3.9863) t + 288.9414 
w _w 

t + 68.12963 
w 


508929.2 


i-c A 

-'1^ ^101325 10330 

|l-(2.11 - 0.053t^) (10"^)j g/cm. 


(45) 


In this expression, t^^ Is the water temperature In °C, C Is the 
compressibility of water (C « 471 x 10“^^ Pa"^), B Is the barometric pressure 
In Pa, I Is the depth of Immersion of the sample In water In mm, and D Is the 
number of days elapsed since the water was boiled. After the measurement In 
water, the specimen Is dried and weighed in air with the same balance. The 
mass density of air is calculated from the following expression given in 
Reference 25: 


r464.56B - H(0.085594t^ - 1.8504t + 34.47) 


t + 273.16 
a 


X 10 


-6 


g/cm^. 


(46) 


where H is the relative humidity In percent and t^ is the air temperature in 
°C. 

Nlchrome (80Ni-20Cr) wire, 0.013 mn In diam, was used to suspend the 
specimens for the weighings In water. The wires were baked at 700-800°C for 
1 h in air at 0.1 Pa (1 x 10“^ Torr) to prepare their surfaces. The oxidized 
surfaces of the wires developed microfissures when the wires were cooled. 
These microfissures cause the water meniscus and surface tension to be nearly 
the same at each point along the wire when it is raised and lowered in the 
water during the weighings. No surface film on the water is tolerable, and 
microscopic dust particles in the meniscus region can invalidate the weight 
measurements in water. 
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Tables 12, 13, and U summarize the results of mass density measurements 
for Ingots L3, L4, and L6. The distance from the tip given in the tables is 
the distance to the midplane of a slice. 


TABLE 12. AVERAGE ALLOY COMPOSITIONS FOR 
VARIOUS POSITIONS ALONG THE AXIS 
OF INGOT L3. 


Sample 

number 

Oiitanoe 
from tip 
(cm) 

Matt 

dentity 

(g/cip3> 

Average 

compotition 

ix) 

L31 

0.25 

7.0319 

0.472 

L32 

0.64 

7.2416 

0.378 

L35 

2.07 

7.5906 

0.221 

L310 

5.01 

7.6312 

0.203 

L312 

6.35 

7.6237 

0.207 

L313 

7.27 

7.6296 

0.204 

L314 

8.22 

7.6285 

0.204 

L321 

12.60 

7.6328 

0.202 

L322 

13.32 

7.6377 

0.200 

L327 

15.91 

7.7200 

0.163 

L328 

16.60 

7.8468 

0.106 


QP0347M'13 


4a 





TABLE 13. AVERAGE ALLOY COMPOSITIONS FOR 
VARIOUS POSITIONS ALONG THE AXIS 
OF INGOT L4. 



SampI* 

number 

Distance 
from tip 
(cm) 

Man 

density 

(p/cm^) 

Average 

composition 

lx) 


L4TIP 

0.29 

7,1098 

0.437 

L41A 

0.73 

7.4947 

0.265 

L41B 

0.81 

7.5886 

0.222 

L41C 

1.26 

7.6051 

0.215 

L41D 

1.53 

7.6144 

0.211 

L41E 

1.79 

7.6220 

0.207 

L41G 

2.33 

7.6174 

0.209 

L41I 

2.86 

7.6176 

0.209 

L41K 

3.39 

7.6194 

0.208 

L42 

4.32 

7.6252 

0.206 

L43 

5.40 

7.6162 

0.210 

L44 

6.43 

7.6211 

0.208 

L45 

7.45 

7.6334 

0.202 

L46 

8.47 

7.6404 

0.199 

L47 

9.49 

7.6327 

0.202 

L48 

10.52 

7.6320 

0.203 

L49 

11.54 

7.6206 

0.208 

L410 

12.32 

7.6154 

0.210 

L411A 

13.07 

7.6171 

0.210 

L411B 

13.33 

7.6263 

0.205 

L411D 

13.87 

7.6474 

0.196 

L411P 

14.40 

7.6976 

0.173 

L411G 

14.67 

7.7458 

0.151 

L411H 

14.93 

7.8229 

0.117 

L411J 

15.21 

8.0653 

0.007 
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TABLE 14. AVERAGE ALLOY COMPOSITIONS FOR 
VARIOUS POSITIONS ALONG THE AXIS 
OF INGOT L6. 



Sample 

number 

Diitance 
from tip 
(cm) 

Matt 

density 

Ig/cm^l 

Average 

composition 

(X/ 


L6A8 

0.80 

7.1005 

0.442 

L6A10 

1.00 

7.1196 

0.433 

L6A12 

1.20 

7.2979 

0.353 

L6A17 

1.65 

7.5192 

0.254 

L6A20 

2.00 

7.560C 

0.235 

L6A24 

2.40 

7.5843 

0.224 

L6B135 

3.50 

7.6069 

0.214 

L6B148 

4.85 

7.6221 

0.207 

L6B162 

6.25 

7.6239 

0.206 

L6C2 

7.75 

7.6067 

0.214 

L6C5 

9.25 

7.6158 

0.210 

L6C7 

10 25 

7.6278 

0.205 

L6DT122 

12.25 

7.6247 

0.206 

L6ET138 

13.80 

7.6249 

0.206 

L6ET147 

14.75 

7.6428 

0.198 

L6ET157 

15.75 

7.6802 

0.181 

L6ET162 

16.20 

7.7101 

0.168 

L6FT169 

16.90 

7.8925 

0.085 

L6FT172 

17.20 

7.8985 

0.082 

L6FT177 

17.70 

8.0522 

0.013 


OM3 07S015 


4.4.2 Infrared Transmission Measurements 

From Ingots L3, L4, and L6, thin ("200-300 gm) slices were prepared for 
evaluation of the radial compositional variation by Infrared transmission-edge 
mapping. The transmission measurements were made at room temperature through 
1-mm dlam areas at regularly spaced locations on a slice. 

The initial measurements on slices from Ingots L3 and L4 showed nc 
Infrared (IR) transmission. Subsequently, the thicknesses of these slices 
were reduced to 20-50 gm by polishing and etching, and the measurements were 
repeated. However, even for these small thicknesses, the slices remained 
essentially opaque to IR radiation. This result was attributed to the highly 
p-type character of the slices. 

The slices from ingot L6 showed up to 40X transmission above the cut-off 
wavelengths. The axial positions along the ingot of the slices on which IR 
measurements were made are shown in Figure 21 (a) . The radial measurement 
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positions are shown in Figure 21(b). The diameter of the circles and their 
positions correspond to the width and relative locations of the incident 
infrared beam. Typical transml jsion data for slice A26 are shown in Figure 
22. For slices A16, A26, E3, E39, and E34, the radial variation of the cut-on 
wavelengths and the corresponding alloy compositions are shown in Figure 23. 
Figure 24 shows the radial variation of the wavelength at 50% transmission. 

The radial compositional profiles indicate concave solid/liquid interfaces for 
the entire growth length and thus significant radial temperature gradients in 
the \/iclnity of the Interface. 


(a) 



Gf>a3-0750'43 


Figure 21. (a) Location of wafers along the growth axis of ingot L6 and (b) the IR measurement 
positions for a wafer. 


( 1 ) 




Wavelength 



Figure 22. Typical transmission spectra for a wafer (A26) from ingot L6. 
Numbers of curves refer to locations on wafer as shown in 
Figure 21. 
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Figure 23. Radial variations of the cut-on wavelength and Cd composition for ingot L6. 
Left-side plots are for areas 1 through 9 and right-side plots are for areas 
10-13, 5, and 14-17 as shown in Figure 21. 
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Figure 24. Radial variation of wavelength at 50% relative 
transmission across slices from ingot L6. Left- 
side plots are for areas 1 through 9 and right- 
side plots are for areas 10-13, 5, and 14-17 
as shown in Figure 21. 


4.5 Alloy Coppositlonal Variation Along the Growth Axis 


4.5.1 Theoretical Modeling 

During directional solidification under conditions of complete mixing and 
near-equilibrium growth, the compositional variation along the growth axis is 
given by 

C (£) = k C (1 - , (A7) 

s o 

where Cg(£) is the solute (CdTe) concentration at a distance i along an ingot 
of length L, is the initial solute concentration in the melt, and k is the 
equilibrium distribution coefficient at the liquid-solid Interface 
temperature. The distribution coefficient is given by k « x^/x , where x® and 
x^ are respectively the solute compositions of the solid and the liquid at the 
Interface temperature. In the derivation of Equation (47), k is assumed to be 
independent of the interface temperature. 
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Under conditions of elemental diffusion In the melt and negligible 
mixing, the compositional variation along the growth axis Is dependent on the 
relative values of the Interface translation rate, R, and the solute diffusion 
rate, D, In the melt. For k > 1 (k < 1), the solute distribution profile Is 
expected to exhibit an Initial transient of progressively decreasing 
(Increasing) solute concentration, a region of constant composition equal to 
the original melt composition, and finally a region of rapidly decreasing 
(Increasing) solute concentration. Smith, Tiller, and Rutter^ have examined 
the mathematical details of the problem and obtained analytical solutions for 
the transient compositional variations in terms of the material and growth 
parameters. Their solution for the initial composition-transient for constant 


R and k 1® given by 



Cq j 

1 


Cg(t) - 1 1 + erf 

X 

2 

(r) J 

. / 2k - 

A, 


X erfc 1 2 

7 

i^) J 


+ (2k - l)e 


-k(l-k) (R£/D) 


(48) 


For the final composition-transient, they obtained the expression 



(2n + 1) 


(1 - k) (2 - k) ... (n - k) 
(1 + k) (2 + k) ... (n + k) 


^ -n(n+l) [R(L-A)/D] 
X e 


(49) 


In the limit Rf » D, the widths of the initial and final transients are 
small, and Cg(l) is essentially independent of I and equal to C^. For D 
and k > 0, the axial compositional profile has an Initial transient such that the 
first-to-freeze solid has a composition Cg(o) * kC^, and with Increasing I 
the composition decreases to a steady-state value of Cg(Jl) “ C^. 

A final transient of rapidly decreasing composition occurs as the melt Is 
consumed. The length of the central steady-state region can be Increased at 
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Iniertace segregation coefficient, k,. 


prior to growth initiation and the uncertainties introduced by the tapered 
geometry of the ampule tips, renders the initial-transient data unsuitable for 
quantitative evaluations. 
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figure 25. The variation of the interface segregation coefficient and interface translation 
rate along the growth axis of ingot L3. 
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Figure 26. The variation of the interface segregation coefficient and interface translation 
rate along the growth axis of ingot L4. 
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Figure 28. A comparison between the measured and calculated compositional profiles 
for ingot L3. The calculations assumed a diffusion coefficient of 
6.0 X 10*5 cm2/s and included the variations in k and R along the growth axis 
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Figure 29. A comparison between the measured and calculated compositional profiles 
for ingot L4. The calculations assumed a diffusion coefficient of 
6.0 X 10'^ cm^/s and included the variations in k and R along the growth axis. 
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Figure 30. A comparison between the measured and calculated compositional profiles 
for ingot L6. The calculations assumed a diffusion coefficient of 
4.5 X 10'^ cm^/s and included the variations in k and R along the growth axis. 
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Figure 31. A comparison between the experimental data for ingot L3 and calculated 
results for various values of the diffusion coefficient. 


ORfG?\’?;3. p ; - . 

O'" POQo 



44 

I 

, 1 ' 

■■ 


ORIGINAL PAGE IS 
OF POOR QUALITY 




Figure 32. A comparison between the experimental data for ingot L4 and calculated 
results for various values of the diffusion coefficient. 
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Figure 33. A comparison betvMtn the expetimental data for ingot L6 and calculated 
results for various values of the diffusion coefficient. 
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Figure 34. A compariion between the experimentel deta for ingot L3 and calculated 
reiult* for a diffusion coefficient of 6.0 x 10"5 cm2/$ and constant k, R. 
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Figure 35. A comparison between the experimental data for ingot L4 and the calculated 
results for a diffusion coefficient of 6.0 x 10'^ cm^/s and constant k, R. 
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Figure 36. A comparison betwraen the experimental data for ingot L6 and the calculated 
results for a diffusion coefficient of 4.5 x 1(T® cm^/s and constant k, R. 




5. THEORETICAL MODELING OF CHARGE-CARRIER CONCENTRATIONS 
AND ELECTRON MOBILITIES 


Theoretical models and coi5*»'in " \jrograms specific to Hgj_j^CdjjTe were 

developed for calculations of carrier concentrations. Hall 

coefficients, Fermi energies, nJ electron mobi’lties as functions of x, 

temperature, and lonized-defect and neutral-defect concentrations. As in 

previous work on the HgCdSe alloy system,^* the Kane three-band model^^ is 

used to describe the Hg2_jjCdj^Te band structure for the energy range of 

interest, i.e., the lowest-lying conduction band (F,) and the two highest- 

o 

lying valence bands (F ) . Higher-band corrections considered in Reference 9 

O 

are neglected because such corrections are estimated to be small for the 
positive band-gap alloys that are of primary concern in this study. The 
electron mobilities are calculated in terms of a microscopic theory of 
electrical conduction derived from the solution of the Boltzmann equation for 
the perturbed steady-state electron distribution function. 

5. 1 Energy-Band Structure and Band Parameters 

The secular equation describing the conduction, light-hole, and spli:- 
off-valence bands Is^^ 

E^ + (A - Eg) E^ - (E A + P^k^) E - -| AP^k" = 0, (54) 


where E is the band energy in terms of the crystal momentum wave vector Jc, A 
is the F^^ spin-orbit splitting. Eg is the F^ - Fg energy b£.id-gap, and P is 
the momentum matrix element between the F^^ valence-band and 
the F^ conduction-band wave functions defined by Kane.^^ The energies are 
referred to the top of the valence bands. The conduction, c, band and the 
light-hole, ilh, band densities of states are given by 


(E) 


k! dk ^ _L 

^2dE 2,2 \ 


(i = c, £h), 


(55) 
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where 
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X (C) 

C 


3[26^C +6 (3 - + -| (1 - 6)5 - I ] [5(5 - D(65 + 1)]^'"^ 

2(|S5 + 1)^^^ 

( 56 ) 




3[2«^5^ - 6(3 - 6)5^ +-| (1 - +t1 ~ 


2(f S? - 1] 


5/2 


(57) 


5 



and 


(58) 


6 - E /A. 

8 

The crystal momentum as a function of energy Is 


(59) 


k 


(I) 


1/2 ^ 


8(5), 


(60) 


where 


8 ( 5 ) 


5(5 - I) (65 + 1) 

(| + 1) J 


1/2 


(61) 


The conduction-band wave functions are given by 

|l[, c, [a |lSa_^±b|(X7 iY) a^)>+ c|Za_)>], 


(62) 



66 


where x. y, and z are the basis set of referred to a coordinate system 
with the x-axls along tc, S is the wave function, are Pauli spin 
functions for spin parallel (+) and anti-parallel (-) to and a, b, and c 
are functions defined by 


[e(6C + 1) (6C + |)]^^^/N, 

(63) 

/2 (C - 1)^^^/3N 

(64) 


and 


c - (C - 1)^^^ («C + |)/N. (65) 

where 

1/2 

N » [C (6C + 1) (6C + |) + I (C - 1) + (^ - 1) (6? • (66) 

The heavy-hole, hh, band is represented by a simple parabolic band given by 


\h 


hiki 

2“o\ 


(67) 


where w is the effective-mass ratio. 

V 


5.2 Calculation of the Temperature Dependence of the Carrier Concentrations 

The Fermi energy, Ep, electron concentration, n^, light-hole 
concentration, and the heavy-hole concentration, n^j^i are calculated from 

the numerical solution of the charge neutrality equation. 


n 

e 


fh 


- n 


hh 




( 68 ) 


where 
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2tt 


n 


Ih 




3/2 




3/2 

1^ 


00 

L 




Q 

/ 


dy ^,h(By)f^(y.-z). 


I 



2p m k_T\ 
V o B \ 

/ 


3/2 


Fi/2(-z), 


(69) 


(70) 


(71) 


z = E^/fc„T 
F B 


(72) 


and 

f^(y,z) = + 1)"^ , (7A) 

In Equation (71), Fermi function of order 1/2, and Np are 

respectively the number of acceptors and donors per unit volume, kg is the 
Boltzmann constant, T is the absolute temperature, m^ is the free electron 

mass, and y = E/k T. 

B 

A literature search was performed to compile the best available band 
parameter values, and these band parameters along with the various material 
parameters required for the mooillty calculations are summarized in Table 
15. Using the values for the band parameters given in Table 15, the 
temperature and x dependences of the Fermi energies and carrier concentrations 
from 4.2 K to 300 K were calculated from Equations (67)-(70) for various net 
donor concentrations. The calculated results agree well with available 
data. A typical comparison between calculated and measured electron 
concentrations is shown in Figure 37 for an alloy composition of x » 0.183 and 
a net donor density Np “ ^A “ * 10^^ cm”^. 


68 



Electron concentration, n (cm' 


3x 10’® 


1 X 10 


16 


3x 10 


15 


1 X 10 


15 


3x 10 


14 


OF POOR Quality 




n^, experimental- 



• ng, theoretical 


No - N^= 5.6 X 10^'* cm'^ 


1 X 10 


14 


0 50 100 150 200 250 300 

Temperature (K) 


Gpojtirso^r 


Figure 37. Theoretical and experimental conductiorvelectron concentration for KIqq g^yCdQ ^ggTe. 



TABLE IS. PARAMETERS FOR ELECTRON-MOBILITY AND CHARGE-CARRIER-CONCENTRATION 
CALCULATION FOR Hg,.^Cd„Te ALLOYS. 


Parameters 


Value 


Reference^ 


Spin-orbit splittmg^®^ 
Heavy-hole mass^^^ 
Deformation potentials^*^^ 


Longitudinal sound velocity^*^) 
Transverse sound velocity**^^ 
Hg mass 

Cd mass 
Te mass 

Lattice constant 
rg-Ps energy gap 


Momentum matrix element 
coupling conduction and 
valence bands!®) 

Transverse effective charges 

Reduced masses 

Dielectric constants!^) 


HgTe LO phonon frequency 
HgTe TO phonon frequency 
CdTe LO phonon frequency 
CdTe TO phonon frequency 


Aq= (-0.18x-H.08) eV 
"'hh^'^0 ^ 

Eq = 11 eV (longitudinal mode) 
g. = 9 5 eV^ 

' V (transverse modes) 

E 2 = 3 eV J 

V|_ = (2.52 0.54x) X 10® cm/s 
vj = (1.60 + 0.27x) X 10® cm/s 
"^Hg = 200.59 a.m.u. 
mQ(j = 1 12.4 a.m.u. 
mjg = 127.60 a.m.u. 
a^ = (0.6460 0.0021 x) nm 
Eg = [-0.03 -*• 5 X 10*^ -r (1.91 - 10'^) x] eV 
(0.135 «:x< 0.203) 

Eg = [-0.303 + 1.73x + 5.6 X 10 (1 - 2x) T 

+ 0.25x'*] eV (0.23 < x < 0.61, 10 K < T < 300 K) 
Eg = [-0.25 X 1.59 X -r 5.233 X lO '* (1 - 2.08x) T 

-I- 0.327x3] eV (0.13 < x < 0.6, 20 K < T < 300 K) 
P= 8.55x lO'SeV cm 


( 1 ) 

( 2 ) 

(3) 

(4) 


(5) 

( 6 ) 

(7) 

( 8 ) 


HgTe: ej/e* = 2.96 
CdTe: ej/e* = 2.35 

'"HgTe = 1-295 X 10-22 g 
'"CdTe = 9-®23 X 1o 23g 

eg = 20.206 - 15.1 53x-^ 6.59091x2 (g) 

-0.951826x3 

£ 00 = 15.1077 - 13.8823x-r 9.88928x2 
-3.67133x3 


wHgTe 

LO 

= (139 -9x) cm'1 


LO 

= (118 -H2x) cm'1 


^CdTe 

LO 

= (154 + 16x) cm'1 

(9) 

cuCdTe 

TO 

= (152 -9x) cm'1 
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TABLE 15. (Concluded) 


(a) Aq = Tg — at k = 0. A linear interpolation was made between the HgTe value and the CdTe value. 

(b) X » 0, 7one center 

(c) The El value is for CdTe (i.e-, x = 1). The variation of Ei with x has not been considered. The value of £3 
IS estimated based on a hydrogemc approximation. 

(d) v|_2 = Cl i/p, v-f^ = C44''P./^ = (~2.3x + 8) g/cm^ (This approximate equation was obtained using data 
from J. Blair and R. Newnham (1011 . The value of ci 1 for HgTe (300 K) Is 50.8 GPa (1 1 ), while the value 
for CdTe is 53.5 GPa (12, 13). The assumption of a linear variation for ci 1 as a function of x yields the 
result. Cl 1 = (2.7x + 50.8) GPa. The value of C44 for HgTe (300 K) is 70.5 GPa (12), while the value for 
CdTe IS 19,9 GPa (12). Assuming a linear variation for c^^ as a function of x, c^^ = (-0.60 x + 20.5) GPa. 
These values for P. ^i 1. and yield the composition dependence for the velocities given in the table. 

(el The value of P was calculated from the equation h^ym^' * (2P^/3) (2/Eg) + 1/(Eg + A)] .where x = 0.204 
and values of m^* = 4,66 x 10*3 rn^. Eg - 61.7 x 10*3 eV, and A= 0.96 eV were used (14) 

(f) These approximate equations were obtained using the referenced data. 

1. A. Moritani, K. Taniquchi, and C Hamaguchi, J. Phys. Soc Jap 34,79(1973). 

2. V.l. IvanovOmskii, B.T. Kolomiets, Yu F. Markov, A. Sh. Mekhtiev, and K.P. Smekalova, Phys. Status 
Solid! 32. K83 (1969), Sov. Phys. Semicond. 1, 1203 (1967). 

3. A. Jedrzejcak and T. Oieff, Phys. Status Solidi B 76. 737 (1976) 

4. D. L. Rode, Phys. Rev. 6 2. 4036 (1970) 

5. J.C. Woolley and B. Ray, J. Phys. Chem. Solids 13, 151 (1960) 

6. J.D. Wiley and R.N. Dexter, Phys. Rev. 181, 1181 (1969) 

7. M.W. Scott, J. AppI Phys. 40. 4077 (1969). 

8 J L. Schmit and E.L. Stelzer, J. AppI Phys. 40, 4865 (1969). 

9. J. Boors and F. Sorger, Solid State Comm. 10, 875 (19721, F Sorger, Dissertation, Universitat Freiburg (1972). 

10. J. Blair and R. Newnham, Metal lurgy of Elemental and Compound Semiconductors (Wiley Interscience, 1967), 
p 399 

11. J.G Mavroides and D F. Kolesar, Solid State Comm. 2, 363 (1964). 

12 D G. Thomas, J. AppI Phys. 32, 2298 (196D, 

13. H.J. McSkimin and D.G. Thomas, J. AppI. Phys. 33, 56 (1961) 

14. M.A. Kinch and D.D Buss. The Physics of Semimetals and Narrow-Gap Semiconductors . Ed. by D.L Carter 
and R.T. Bate, (Perg ion Press, 1971), p. 461. 

t A good compilation ot data on the properties of Hg.j ^Cd^Te is given by R Dornhaus and G. Nintz, Springer 
Tracts in Modern Physics, Solid State Physics. (Springer Verlag, 1976) Vol 78. p. 1 ff. 
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5.3 Theoretical Modeling of the Electron ^foblllty 

The electron nobility calculation Includes the following scattering 
mechanisms: longltudlnal-optlcal phonon (LO), longitudinal- and transverse- 

acoustical phonon (ac), heavy-hole (hh), and alloy disorder potential 
(dls.). The extrinsic scattering mechanisms Include charged (11) and neutral 
point-defects (nd). The current densities, conductivities, and mobilities are 
calculated from the perturbed electron distribution function f(t) given by 

f(^) = f - kc'(E) f' cos\(i, (75) 

o o 

df 

where f„ is the unperturbed distribution function and f' “ quantity 

is the angle between the wavevector tc and the applied electric field. The 
pertubation function c' (E) depends only on the energy. The current density is 
given by 

00 

J = / l^^c'(E) f'dE, (76) 

^ 37t2^ J 


where e is the electronic charge, -fi is Planck's constant, and c'(E) is given 
by the solution of the Boltzmann equation. 


3f 


!l* ■?(*), •(§;), <■ ■ “■ 


ac, eh, dis., il, nd) . 


(77) 


The sum 


■? (tt). • 


is the rate of change of the distribution function caused 


by the various scattering mechanisms, and^'lYj is the rate of change of the 
distribution function caused by the application of a static electric field. 
For the steady-state. 


Miii 

3t 


- 0 


(78) 
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The field term is given by 

'cost, 


/3f\ d" e , 

(hL ■ ^ ^o‘ 


(80) 


where S“ is the applied electric field. The scattering term is given by^® 


2(17) • -7^ S tkc'(E) -k'c'(E')] dk'costi;, 


(81) 


where 


V. (k, k') = W. (k', k) f (k') [1 - f (k)] - V, (k', k). 
1 10 01 


(82) 


W^(k', k) is the transition probability per unit time per unit volume and is 
given by 


(k',k) 


H^(k', k) 




(83) 


where H^(k', k) is the matrix element for scattering from state k to k' for a 
given scattering process, is a delta function, and Eg is the energy absorbed 
or emitted in the scattering process. 


Equation (79) reduces to a linear finite-difference equation in c'(E), 
and the conductivity is determined from Equation (76) by using the variational 
method of Kohler^® as modified by Howarth and Sondheimer^^ and Ehrenreich.^^ 
The pertubatlon function is expanded in a complete set of trial functions. 



c 4 . (E), 
n n 


(84) 
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and the are determined by the requirement that c'(E) be a stationary 

28 

point of a certain conserved Integral. 

The various scattering mechanisms and their contributions to the 
Boltzmann equation are summarized in Appenoix A* The computer program for 
calculating the Fermi energy, the charge-carrier concentrations, and the 
electron mobilities is given in Appendix B. 

5.4 Calculations of the Temperature Dependence of the Electron Mobility 

Calculations were performed to ascertain the relative contributions of 
the various scattering mechanisms to the conduction electron mobility as 
functions of temperature from 4.2 K to 300 K for alloy compositions near 
X ■ 0.2. The material parameters used in the calculations are given in Table 
15. Figures 38 and 39 show the calculated results for the sample described In 
Figure 37. The curves designated by the symbols » ^eh’ ^11 

correspond to the calculated mobilities for the individual scattering 
mechanisms. Tlte calculations shown in Figure 38 were made assuming that all 
ionized defects are donors, l.e., the donor density is equal to the residual 
conduction-electron density at low temperatures (T ^ 50 K) . The calculated 
results in Figure 39 assumed a donor density of Np * 7.6 X 10^^ cm”^ and an 
acceptor density of 2.0 X 10^^ cm”^. The calculated results indicate that for 
alloy compositions of x = C.2, longitudinal optical-phonon scattering is the 
dominant mobility-limiting mechanism at high temperatures 75 K), while 
lonized-defect scattering is the dominant mechanism at low 
temperatures CS 75 K). 

The agreement between the temperature dependence of the calculated total 
mobility ,ii^, and the experimental data is reasonably good. However, as the 
temperature Increases, the calculated and measured mobilities systematically 
diverge until at 300 K the calculated mobilities are about a factor of two 
higher than experimental values. A probable reason for this discrepancy is an 
overestimate of the screening of the electron-LO-phonon interaction. The 
approximate screening function used in the ilculations treats the electrons 
as always screening and Ignores the possible antiscreening effects considered 
by Ehrenreich^® for small momentum transfers. The calculations also neglect 
possible scatterings by neutral point-defects and by line-defects. These, 
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combined with tho uncertainties in the various material parameters, could 
explain the differences between the calculated and measured results. 
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Figure 38. A cotnpariion of the temperetiire dependence of the calculated and measured conduction-electron 

mobilities for an alloy composition of x > 0.183. The calculated results are for a donor concentration 
equal to the residual election concentration in Figure 37. 


75 



ORIGINAL PAGE 1 % 
OF POOR QUALITf 



C- 

Figure 39. A comparison of the temperature dependence of the calculated and measured conductioirelectron 
mobilities for an alloy composition of r = 0.183. The calculations assumed a donor concentration 
of 7.6 X cnr3 and an acceptor cor.rentration of 2 x 10^^ cnr3. 
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APPENDIX A; ELECTRON SCATTERING MECHANISMS 
AND THE SOLUTION TO THE BOLTZMANN EQUATION. 


The mathematical expressions used for the various scattering mechanisms 
are essentially the same as those used for the calculations of electron 
mobilities in the HgCdSe alloys described In References 9 and 10. Some of the 
pertinent mathematical details are included as part of this report to aid in 
interpretation of the computer programs. 

A. 1 Longitudinal-Optical (LO) Phonon Scattering Mechanism 

The strength of the LO-phonon/electron interaction is proportional to the 

•k 

Callen effective charge, e , associated with relative displacements of atoms 

c 

in a unit cell. The Callen effective charge is related to the Szigeti 
k k 

effective charge, e , and the total transverse charge, e , by the relations 
o X 


* r ' 

i s - 

C 


+ 2 


3 e 


(Al) 


and 



(A2) 


where is the high-frequency dielectric constant. The interaction is 
assumed to be screened by the dielectric function given by 

(k® 

q) - , (A3) 

q 

where is the longitudinal-phonon frequency, q is the phonon wave vector, 
e £h 

and ky^ and ky^ are the electron and light-hole Fenni-Thomas momenta given by 
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,,«h,2 2e 
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^^2 / 2 „-„ . , y /2 - 
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(A5) 


where 

p' = 3h^ E /4m P^, (A6) 

g o 

and the other parameters have been defined in Sections 5.1 and 5.2. 

The scattering tern of the Boltzmann equation for the LO-phonon-electron 
interaction is given by 



A-n^k^KF, . 

2 2 
"h ( 1 ), X s E e 
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(c') f' cost|», 
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where 
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M, 


(i 


HgTe, 1=2- CdTe), 


(A8) 


is the number of unit cells of the i^*' type per unit volune, is the 
reduced mass of the i^^ unit cell, and 


iJ? (c') 
LO 


»■* <"p * 

1 

+ (f /f ) n (c' R •• c'S ) I . 
-op-- - 


(A9) 


The subscripts (±) refer to the evaluation of the function at y ± 0, where 
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y - E/k„T and 0 = -ftu, /k„T. 
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The function, np, is the phonon occupation number given by 




n = (e - 1) 
P 


(AlO) 


The functions and are defined as follows: 
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s^(y, 0) = ^ ^n ^y* 
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* [(6By +|)/N(6y)] [(| 6By + 1)/By(66y + 1)]^^^, 
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P* (y, 0) = (aa^)^ - bb^[bc^+ cb^) //2] 


-by (bc^ -b cb^)^ +j (bb^)^. 
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P* (y. 9) “ 2aa^ 
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(A24) 


P^. = 


= ss^ . 


and 


2 -B ,,,e .2 , ,.Jlh-2, 

“ 2u'm E e ^ 

o g oo 


(A25) 


(A26) 


The s)rinbols a, b, and c have been defined In Section 5.I. In the above, 

R and S are zero for y-9-6 ^<0. 

A. 2 Acoustical Phonon (ac) Scattering Mechanism 

The longitudinal and transverse acoustical-phonon/electron interactions are 
characterized by three deformation potentials: for the longitudinal mode 

and Ej and E 2 for the transverse modes. The acoustical-phonon contribution to 
the Boltzmann equation is given by 



2k„T E - 

B o , , , 2 _ . 
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irh p V, 

m L 


(A27) 


X 



cost c' (E) , 


where Pjj^ls the mass density of the material, and Vj and v.^ are the 
longitudinal and transverse velocities of sound, respectively. 

The longitudinal and transverse scattering functions, Fl and F,j,, are given by 


A? + I (AjA2 + A^, 

5^ ^2^3 "*"3 ^3 


I (2AjA3 + A 2 - A^) 


(A28) 
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where 


Aj = a^ + (Ej/Eq)(| b^) + (E 2 /Eq)(c^ + | b^), 

^2 “ ■ (E^/EqHc^ +1 b^) + (E^/EpXSc^ +1 b^), 

A 3 = 4(E2/Eq)(c2 -\hh. 

R - i 

»! - 4 b , 

B2 = -b^ [b(b + /2 c) + i b^ - (b + c//2)2], 

= b^[(b + /2 c)^ + b(b + /2 c) - 2(b + c//I) (b + /2 c) 
- (b + c//2)^] , 

B^ •= 2b^ (b + c//2) (b + /2 c) , 
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1 f, 2 2 . 1 , 2 1 ,2v2, 

-j [b c + -j (c - -2 b ) ] . 


(A39) 


A. 3 Ionized Impurity (11) Scattering Mechanism 

The electron/lonlzed-defect Interaction Is assumed to be screened by the 
dielectric function 


e(q) - 

q 


where k^ll Is the heavy-hole Ferml-Thomas momentum, given by 
r 1 


(khh^2 

^ FT 


ei 


3/2 


"- 1/2 


(A40) 


(A41) 


Is the low-frequency dielectric constant, and ^^ 1/2 Fermi function 

of order -1/2. 

The contribution to the Boltzmann equation from lonlzed-lmpurlty scattering Is 



.(NpZg . M^Z^) 

h E s^ 
o g 


>11 f' cosmic' (E) , 


(A42) 


where and are the donor and acceptor concentrations, respectively, and 
Zp and are their charges in electron-charge units. The scattering 
function, is given by 
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A. 4 Electron-Hole (eh) Scattering Mechanism 

The contribution to the Boltzmann equation from electron-hole scattering 
is 
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TT n, . e X 
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Ti E E s 
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where is the same as for ionized defect scattering except that g^ is 
replaced by l.e.-. 


j^2.,,e v2 . ,,ih.2 . ,,hh»2, 

2y'm E e 

O g 00 
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A. 5 Compositional-Disorder (dis») Scattering Mechanism 

Electron scattering by the random short-range potentials resulting from 
alloy compositional disorder is treated as scattering by a random distribution 
of square-well scattering centers, which have unit-cell dimensions and depths 
approxlinately equal to the difference in band gaps between CdTe and HgTe, The 
contribution to the Boltzmann equation from this process is given by 



2 * 
— §■ (ii'n ) sX — 


(1 - x) E 


dls. 
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♦ c'f'costfi, 
dls. o 


(A53) 


where is the number of unit cells per unit volume, Is the disorder 
potential, approximately equal to the difference in CdTe and HgTe band gaps, 
and is a wave-function-overlap Integral given by 
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A. 6 Neutral-Defect (nd) Scattering Mechanism 

The scattering of electrons by neutral point-defects Is also described in 

terms of scattering by a random distribution of square-well potentials. Thus 

the contribution of neutral point-defects to the Boltzmann equation is also 

2 2 2 

given by Equation (A53) with x(l-x) E,. /N replaced by N w E , where 

dis. a n n n n 

is the density of neutral point-defects,i»jj is the volume of the potential 
well, and E^ is the strength. 

A. 7 Calculation of the Electron Mobility 

The total scattering contribution to the Boltzmann equation from the 
scattering mechanisms considered is 
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with the following defined quantities: 
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The substitution of Equation (A55) ^nt, ini; st idy-state Boltrmann equation. 
Equation (79), yields the following finite-difference equation: 
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To obtain a solution, the perturbation function, c'(E), expanded In a 
complete set of trial functions. 
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(A63) 
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'0 


and the r.^ are determined by the requirement that c'(E) be a stationary 
point of a certain consen'_J Integral. The set. 


(y - rh" . 


(A64) 


is used for the trial functions, and the expansion is carried to third order 
(i.e., n * 2). The current density Equation (75), is given by 
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The electrical conductivity, o. Is given by the relation 


o - J^/<? , (A70) 

and the electron siobillty Is given by 

U « o/n e. 
e 


(A71) 


APPENDIX B: COMPUTER PROGRAMS FOR CALCULATION OF CHARGE CARRIER 

CONCENTRATIONS AND ELECTRON MOBILITIES 


B.l Description of Symbols - M a in Program 


Line no. 

Comments 

00330 

Input parameters; TMIN = minimum temperature (K), TMAX = maximum temperature <K), 
DELT = temperature increment, DONOR DENSITY (cm’^), ACCEPTOR DENSITY (cm'3), 
PROPORTION OF CdTe IN ALLOY (at.%), DISORDER ENERGY IN (eV). 

00500 

ES = electron charge in c.g.s. units. 

00510 

HB = Planck's constant (erg s). 

00520 

K = Boltzmann's constant (erg/K). 

00530 

MO = free electron mass (g). 

00540 

PI =7T. 

00550 

CEVE = units conversion factor (eV to ergs). 

00560 

FFCTS is an external program for evaluating Fermi functions. A separate data file of the values 
of Fermi functions is required. 

00705 

X2 = at.% of CdTe in alloy, x. 

00706 

XI = at.% of HgTe in alloy, 1 — x. 

00741 

HCSPOL is an external program used for determining a number of parameters used in the calculations. 

00790 

TMIN = minimum temperature (K). 

00800 

TMAX = maximum temperature (K). 

00810 

DTEMP = temperature increment. 

00820 

ND = donor density. 

00840 

NA = acceptor density. 

00880 

DEV = spin-orbit splitting (in eV) as a function of x. 

00885 

D = spin-orbit splitting (in ergs) as a function of x. 

00890 

UV = heavy hole mass (in electron mass units) 

00960 to 
00980 

EOEV, E1EV, and E2EV are deformation potentials in eV. 

00990 

VL = longitudinal sound velocity as a function of x. 

01000 

VT = transverse sound velocity as a function of x. 

01030 

MOO = Te mass in amu 

01040 

Ml = Hg mass in amu 

01050 

M2 = Cd mass in amu 

01060 

PHI = mass density 

01070 

AL = lattice constant as a function of x. 

OHIO 

EDISEV = disorder energy in eV. 

01120 

EDIS = disorder energy in ergs. 

01140 to 
01160 

EPHIO, EPHM, EPHI2 are deformation potentials in ergs. 


aPOMTSMI 
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Line no. 

Comments 

01170 

E1S = optical phonon scattering contribution. 

01180 

E2S = optical phorwn scattering contribution. 

01190 

M1B = reduced mass. 

01200 

M2B = reduced mass. 

01210 

VA = the lattice constant cubed divided by 4. 

01220 

WL1 = longitudinal optical (LO) phonon frequency (HgTe), kg where = LO phonon 

temperature. 

01230 

WL2 = longitudinal optical phonon frequency (CdTe), kg where Ofj = LO phonon temperature. 

01280 

EG = Eg (ergs) as a function of x and T (temperature). 

01285 

P = momentum matrix element coupling conduction and valence bands (eV-cm) as a function of x. 

01290 

DELTA = Eg divided by the spin-orbit splitting. 

01330 

SECULAR = an external program for finding the limit E = E (k) at the zone boundary. 

01370 

3/2 1/2 

2e2 /2>i'mo\ /'<bT\ 

CON = I 1 1 , a factor that enters into Fermi-Thomas momenta expressions, 

" \ t,2 / \ 0 / 

i.e., (k®-p)^.(k^'})^,and(k^?J-) . 

01780 

SIMP = an external program for evaluating integrals. 

01790 

KFTE2 = electron Fermi-Thomas momentum. 

01810 

KFTLH2 = light-hole Fermi-Thomas momentum. 

01820 

KFTHH2 = heavy-hole Fermi-Thomas momentum. 

01840 

THETA1 =dg^/T (see line no. 01220). 

01850 

THETA2 = /T (see line no. 01 230). 

01860 

GINF2 = g^ a parameter involving Fermi-Thomas momenta. 

01870 

G02 = gg, a parameter involving Fermi-Thomas momenta. 

01890 

A1 = longitudinal optical (LO) phonon scattering contribution. 

01910 

A2 = acoustic phonon (ac) scattering contribution. 

01930 

A3 = ionized impurity (ii) scattering contribution. 

01950 

A4 = electron-hole (eh) scattering contribution. 

01960 

A5 = disorder scattering (dis.) contribution. 

01980 

NP1 = phonon-occupation number. 

01990 

NP2 phonon occupation number. 

02150 

Q1 =i^i (matrix element). 

02160 

Q2 = ;?’2 (ftiatrix element). 


OPOM750-M 
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Lina no. 

Comments 

02170 

03 (nutrix element). 

02200 

QL1 = (SB, (matrix element). 

02210 

QL2 = (SB, 6)2 (matrix element). 

02260 

SIGMA » conductivity. 

02270 

UH = Hall mobility. 

02280 

EGOPT = kgr/BETA. 


GPO3-07S0-«$ 
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B.2 Description of Symbols - Subroutine Programs 


ORIQffJAL PAGE S3 
OF POOR QUALIlY 


Line no. 

Comments 

A. 

SECULAR - The program used to find the roots (E^, E 2 , E 3 ) of the cubic secular equation describing 
the conduction band, light-hole band, and split-off valence band. 

02520 

CUBE - A subroutine for finding the roots of a cubic equation + C 2 *^ + cix + Cq = 0. 

B. 

FUNC - The program used to find the carrier concentration (electrons, light-holes, and heavy-holes) 
as a function of temperature. 

02700 

SIMP - A subroutine used for evaluating a certain class of integrals. 

02710 

NE = electron concentration. 

02780 

NLH = light-hole concentration. 

02800 

NHH = heavy-hole concentration. 

C. 

FV - The program used to obtain the Fermi distribution function 

D. 

FUNC1 - The program used to obtain the product FV and LAMDA. 

03000 

LAMDA - A function used to calculate the conduction band density-of-states. 

E. 

FUNC2 - The program used to obtain the product of FV and LAMDAV. 

03080 

LAMDAV - A function used to calculate the light-hole density-of-states. 

F. 

FVP - This program evaluates the function exp (y + z)/[exp (y + z) + 1 ] ^. 

G. 

FUNC3 • This program evaluates the function \ exp (y - z)/[exp (y - z) + 1 ] 2. 

H. 

FUNC4 - This program evaluates the function Xjj, exp (y - z)/[exp (y - z) 1] . 

1 . 

LAMDA = an expression used in calculating the conduction-band density-of-states. 

J. 

LAMDAV = an expression used in calculating the light-hole density-of-states. 

K. 

S - This program evaluates a function used in determining the crystal momentum. 

L. 

LI - This program evaluates L^”^ 

N. 

NABC - This program evaluates the function Fj^^ (y), which appears in the LO-phonon scattering 
operator 

0 . 

RSL - This program evaluates functions appearing in the LO-phonon scattering operator. 

04390 

ROP1 = pq *' (y,0i> 

04410 

ROP 2 = po'^ (y. 02 ) 

04440 

ROM1 =Po~(v.0i) 

04470 

ROM2 = Po~(y> ®2^ 

04530 

R1P1 =pi'^(y,0i) 

04540 

R2P1 =P 2 ^ (y,0i) 

04590 

R 1 P 2 = p,+ (y,e 2 l 

04600 

R 2 P 2 = P 2 ^ (y, 02 ) 

04660 

R1M1 =pf (y,0i) 

04670 

R2M1 =P 2 "(y,fll) 

04730 

R1M2 = pf (y, 62 ) 


QP0347MHM 
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Comments 


Line no. 


04740 

R2M2 =/;)2~<y'^2> 

04760 

SY2 = $2 

04770 

Q1P = qi'^ 

04780 

Q2P = 

04800 

Q1M = qf 

04820 

Q2M = q2 

04840 

P1P= P/ 

04850 

P2P = P2^ 

04870 

P1M = Pf 

04890 

P2M = P 2 ' 

04910 

V1P = r/ 

04920 

V2P = r 2 ~^ 

04940 

VIM = tf 

04960 

V2M = r2' 

05000 

V0P1 = Vq^ 

05010 

V1P1 = Vi'*' (y,«i) 

05020 

V2F, = V 2 '^ (v.f*i) 

05040 

V3P1 = V3'*^ (y,0i) 

05090 

V0P2 = (y, O 2 ) 

05100 

V1P2= Vi^ (y, 02> 

05110 

V2P2 = V2'^ (V, 02> 

05130 

V3P3 = V 3 '^ (y,02> 

05190 

V0M1 = Vo"(y,0i) 

05200 

V1M1 - Vf (v,0i) 

05210 

V2M1 = V2‘(y, 0i) 

05230 

V3M1 = V3'(y,0i) 

05290 

VOM2 = Vo~(v,02) 

05300 

V1M2 = Vf (y,02) 

05310 

V2M2 = V2‘(y,02> 

05330 

V3M2 = V3'(y,02) 

05390 

RP1 = R+ (y,0i) 

05400 

SP1 =S+ (v.flil 

05430 

RP2 = R+ (y, 

05440 

SP2 = S+ (y, 62 ) 

05480 

RM1 = R_(y, 
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Line No. 

Commentt 

05490 

SMI =S_(v,0^) 

05560 

RM2= R.(y, flj) 

05570 

SM2 = S.(v,» 2 > 

05620 

PIO = ttq (v) 

05640 

PI 2 = 7T2 (v) 

05690 

PHIO = «o (V) 

05710 

PHI1 =0, (y) 

05730 

phi 2 = 02 (y) 

05770 

PHI =4> 

05810 

> 

II 

_> 

05820 

A2 = A 2 (y) 

05830 

A3 = A 3 (y) 

05840 

A4 = A 4 (y) 

05850 

CD 

II 

OD 

05880 

82 = 82 (y) 

05890 

83 = 83 (y) 

05900 

84 = 84 (y) 

05910 

85 = 85 (y) 

05930 

86 = 8 e (y) 

05940 

FLY=FL(y) 

05960 

FTY = Fj (y) 

1 3010 

L 8 = Lb (y) 

P. 

ROOT • This program calculates limits. 

Q. 

INTEG - This program evaluates integrals. 

R. 

FMN • This program calculates the elements of the determinant F^p^ (y). 

S. 

FM - This program evaluates a function involving FVP and LI. 

T. 

HCSPOL - This program calculates as a function of x, the TO and LO phonon temperature, the 
transverse effective charges, the Callen effective charges, and the dielectric constants Cq and Coo- 

07520 

ESTR1 = HgTe transverse charge. 

07530 

ESTR2 = CdTe transverse charge. 

07531 

XM1 = HgTe reduced mass in g. 

07532 

XM2 * CdTe reduced mass in g. 

07533 

A = lattice constant as a function of x. 

07540 

EPS1 = Coo as a function of x. 

07550 

EPSO = 6 q as a function of x. 


QP03-07S048 
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Lift* No. 

Commonti 

07560 

WT2 = w-p (CdTe) as a function of x. 

07570 

WT1 = coj (HgTe) as a function of x. 

07630 

THT1 = HgTe TO • phonon temperature. 

07640 

THT2 = CdTe TO - phonon temperature. 

07^<31 

THL1 - HgTe LO • phonon temperature. 

07782 

THL2 = CdTe LO ■ phonon temperature 

07860 

EC1 = Callen effective charge for HgTe. 

07870 

£C2 = Calien effective charge for CdTe. 


QP034)7S0^ 
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CNF. CN PB 2 

ro Iff - IQ 
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ORiGfNAL PAGu 5S 
OF POOR QUALITY 


Itl 


UV:0.53 

E0EV:9.5 

E1EV=4. 

E2EV-3. . 

V'L:! 2,52*0.E.4»«2)«1 -OES 
VTri 1 .60*0.'.’7»X2 )k 1 -OE5 
HOOr i 27. 60 

, XHI 

Ml =200.59 
M2rl U.i 

^ 

EHI =+.-(X2*112.T"<1 *200. 594 r-'7. 60 )«1. 66056551 !$. 460*0 • 02 1 «X2 

ki <3 

qL = 16.460*.021*S2l*l .OE-8 

E0ISEV=PflRHf7i 

EDI5:CEVE*EDISE\' 

EPH10=rEVE»EOEV 

epmii=:eve*eiev 

EPHI2-CEVE*E2EV 
El S:F1RS»ES 

E2Src 2PS»ES 

^ -I 

llie=l .660'.3E-24»M00«M1 /(MOO*M1 • 

M2B=1 . 66043E-2*KnOO»M2/lMOO*f12 1 
Vfl=BL»BI-«9L24 • 

WL1=K»VHL1/HB 

XNBri2y*PHI /f f100*X l«Ml ♦X2»M2 1 / 1 . 66043E -24 
T = :niN — 




fTG=EEJE«(-0.25*l.S9.X2*0 327«X2«X2.X2*5 . 233E- 4- 1 .0- 2 . 08. :;2 1 *T ) 

P=8. 6E-8«CEVE 
DELTfl=EG/D 
BETB=KT/Ee 
61=1 ./BETA 

UP=3-»«MB»HR.EG/( 4..MQ.P.P' — 

— I 

1 CALL SECUlliRl 2 .»P1 /W. .VBC. VBV.ESI j 


> , 

CONTIN jE I 


rBC=rB0/KT 

rBV:-r{iV/KT I 

C0N=2 . »iUP»M0»<KT/( HB.HBl | 
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I 2=Z*0ELn 

n'OU-FUNCf Z ) I 
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C KKKalXKaiKSKKaiaiKMKKXKMIiKIIMXaOiadIMliMilaiKKKMKKMIlKKMXiIMKKXlIKKKII 


SUBROUTINE SECULfiRI KS .El . E2. E31 


C XKXXXXXXftXXXXXXKKXXXXXXXXXXXXXXXXXXXKXXXXXXXXXXXXXXXXXXX 


REAL KT. NO. MO. NR 
REAL KS 

COMMON /CONSTl/ PI .P. HB.KT .NO . UV. MO.D . EG. YBC. T BV. NR . Z I . ZA 

COrPxPxKSxKS 

C2=0-EG 

C1=-(ECxD4-C0) 

CO=-2.xDxCO/3. 


CALL CUBEIC2.C1.C0.E1.E2.E3) 


1 RETURN 1 
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original 

OF POOR QUALITY 


C MKKKiiKKMKMKKKMKKWMKKKKKKifIKKKXKIk'KXXXkKXKKXKKKKK 

FUNCTION FUNC( X) | 

^ ■ 

C aOIKXKIiaiKKKKXKKKKKKKKXXKKKKKXXXKKKXKKMIlKXKKIiXItK 

T „ - : ':rzL_„ 

REAL NE.NLH.NHH 
REAL KT.ND.NO.NA 
REAL NDD.NDl 

COMMON /CONSri/ PI .R, HB.KT .NO . UV. MO.O . ES, TBC. YBV. MA.ZI ,ZA 

COMMON /NNN/ NE.NLH.NHH.NDl 

COMMON /VAR/ Z.CON 

COMMON /C0NST2/ BETA. DELTA ,BI 

EXTERNAL FUNCl .FUNC2 

^ — ■ 

NE=0.0 
X2 = 10.4-BI 


IF( Z.GT.BI 

i 

X2=10.4-Z 



X2=AMIN1 (TBC.X2) 


T 

i ^ 

GO TO 100^ > 

a 

2 

CALL SIMPI FUNCl .BI , X2, 0 .01 . 10. SI 1 . NE.N . lER! 


^ 

rz~ 

1 NEsCONkNE 


id 


>■ 



rNLHrO.O I 
V 

CENT. CNPG 2 


> 
























I FUNCTION FV( r.Z)~! 
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C »<)(XK><>SKK)ir K»;KK>.KKK»i9(Kl!>;»»KKKKKX>tKKl(X)(KK)(KXKKX 
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3: 


- «K KM K« KK KM «R K« KK X 


RERL FUNtriON LRnOR(yi“ 

“5T^ 


^RRRRRXXRRXRXRXItXRXRXRRKXRRRXXXXXXRRXXRXXXRKRRRRirXRXRXRXR'r 


COMMON /C0NST2/ BETR .OCLTfl, BI 
Br=BETP«Y 



X=3.RDELTflRBr/2. . 

I LRMDRsS. «LRMOfl KSQRr( Sy«( BY-l ■ ) k( DELTRrBV'I-1 . 1 /( XrXkXk XrX) 1/2. 

^ 

RETURN 


Feno 1 


IC »*R* R»<R»«<R*«mH<»(»(RRRXRXK»iKXIIK»(l<KRRmtRRl(.ni«»,MK;HK*>(Klf(XRm(l( 1 

J- 


RERL FUNCTION LRMDRVfT) 


C ***»»»*»»«»«RR»R»(RRRRR»<RI>KRIi'KifX»RliK>iKRXKxR»<RKRm<XI<XI( 




I] 


COMMON /C0NST2/ BETR . DELTA. BI 
BTzBETflKY 


LRMDRV=2.RDELTR«DELTfl»6'«<Br«6y-0ELrfl«(3.-DELTfl)*BrKBY 

»4.rM . -DELTQ)RBr/3.4-2./3. 


Xr3.*DELTR*BY/2.-l . — — 

LRMDRy=3.RLRMDRVi(S0RTtBrR(BrM . )R(DELTflwBr-l. )/C XkXkXkXkX )] /2 . 

— j 


I RETURN I 

IW1 
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*) 




~cND 


ES-X-Qlfl, 
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(4)^ 



I END I 
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CCNT. CNPG 


\t 

















CCNT. CWPG 


10 
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ro imr n 
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C XMKKXXKKKKXKKKIIXKKXKKKKKKXKKKKKalXKKKKKKXKKKaU. 
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SUBROUTINE HCSPOLl X.TMTl . THT2 , 
EPSI.ECl .EC21 


THLl .ThL2.ESTPl ,ESTR2,EPSa. 


ESTR.z2.96 

ESTR2r2.35 

Xm = . .2950" - E-22 

XH2z9. 92339. £-23 

fi:(6. 4604-0 2UX) »1 -OE-3 

EPSIzl?. 10’’7-1 3. 8823*..'', 49 .88928»X«X-3.E71 33KX)KX)i<X 
EPS0:20. 286-15 .1 531 *X*8. 5909 UXkX- 0. 95 1 3 26xX*Xi<X 
WT2=.52.0-9.0*X 


WTlzl 18. 0412. Ox) 

ESz4. 80325E-1 0 
HB=1 .0545S19E-27 
BKzl . 380622E-16 
PI=3. U1593 
Cz2.9979250E10 
THTlz2.»iPI*C*HB/BK)KWTl 
THT2z2. »<PIi«CkHB/8K«WT2 


OnTlzBKxTHTl /HB 
OnT2=BK»TriT2/HB 
F.r4.*PI/(n»*3./4.. 
F2=4.xPI/(fi»<*3./4 . 
S:=F1/0MT1**2. 
S2=F2^0MT2»»2 . 


IxESTRl **2. /XMlxl 1 .->)*ES**2. 
1XESTR2XX2. /XM2*X»ES»»2 . 


THLf12zl ./2./EPSI *( EPSUl THTl **2. 4ThT2».».2 ■ )4S1 »ThT 1»(»2. 
4S2«THT2»(i«2.-SQRr( (EPSIk( THT2*<«<2 . -TMTl **2. )4S2«THT2»*2. 
-Sl«THTl»«2. )*x2-.44.»(THTl*x2.xTHr2«»2.*Sl«S2) 1 


THLP2:1 ./2./EPSIx(EPSIx( THTl **2- 4THT2**2 . )4S1 »THTU»2 . 
4S2»rHT2»x2.4 5QRTl ( EPSI*( THT2**2 . -TMTl **2. l4S 1»THT2»»2 . 
-Six THTl XX2. )x«2.4i.KTHTlxK2.»<THr2»x2.»<Sl»S2) 1 


THLOzSQRTt t EPSI 4S 1 1 /EPS I 

THL12=THLt12 

THL22zTHLP2 


xTH'l 



TO 10 > -lfc 




THL1 2=THLP2 


CENT. CM FG 2 


Efi-LCE Z 
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(t 


THL22=THL(12~j 


0 

1 

m 

Q S 


1 iMLlzSQItT{ THLl 2) 
r 


THL2=SQRr( rHL22) 
rHLl=( THL14-THL2J/2. 
THL2=THL1 
EClzESTRl/EPSI 
EC2=ESTR2/EF Sr 

I ^ 

RETURN 


I END I 
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OF POOR QUALlPf 


OOlOOSnXED 

0(/110$SOURCE»CUtlE«M0 3268 

00120SSOUHCE«EFCTS.M0326B 

00130$SOU«CE«SIMP,M0326A 


1150 


OC 


00 

00 . . 
00160 
00170 

ooieo 

00190 

00191 
00200 
00210 
00220 
00230 
00240 
00250 

0026Q 

00270 

00280 

00290 

00300 

00310 

003?0 

0O33O 

00340 

00350 

00360 

00500 

00510 

00520 

00530 

00540 

00650 

00560 

00570 

00580 

00590 

0U60O 

0U610 

00620 

00630 

00640 

00650 

00660 

00670 

00680 

0069c 

00700 

OU705 

00706 

00710 

0075o 

00730 

00740 

00741 
0074? 

00743 

00744 

00745 

00746 

00747 

00750 

00751 
00760 
00770 
00780 
00790 
ooeoo 
00610 
00820 
00630 
00840 
00850 
00880 
00885 
00890 


MERGE FfCTSi INTRPLt CUBE » SIMP 
PEAL KFTE2iKFTLH2,KrTHH2 
PE AL MN(3«3) tMLl(3» 
logical GENO 

real NE«^LHtNHH*NPl,NP2 
REAL KiNOfKT,MO*M. fM2»M00tNA 
real NUl 
REAL MIB*M28 

DIMENSION ISBSCP<71»VALUE17) »PARM(7) 

COMMON /NNN/ NE»NLH»NHH»N01 

/CONSTl/ PItPtH8fKT»NO»UV»MO»OtEG»YBC»VBV»NA»2I»2A 

/VAR/ ZfCON 
/NM/ NtM 

/C0NST2/ BETAtOELTAfBI 

I,TmEU2»< 


COMMON 
COMMON 
COMMON 
COMMON 
COMMON 
I VL*VT,A2f A3,A4,A5 


/C0NST3/ THETAItTHEU2*GINF2tG02»EPHI0*EPHIltEPHI2t 


'COMMON /COr-ST4/ NPl tNP2fXltX2*AI 
EXTERNAL FUNC3tFUNC4 
external FMN.FM 
W«ITE(6f 1) 


I FORMAT ('Parameter LIST*/* , 
6'OELT*/* 4 •* DONOR OENSII 

— corr IN 


•• 


•• 


V* ^5*»»^ACc|pT0R 


TMAX«/» 3 
OENSITY'7 


•• i< 


ALLOY»/» 7 •'* DISORDER 


OF PARMS TO BE CHANGED AND 

i/ipr 

!CP(1>( 


•/•FOR EACH'CHANGE'//) ■ 
TfVALUEd) #1»1»NCHN6E) 


6' 6 *• PROPORTION Of 

6ENEH0r IN EV'//) 

ES = 4. ,)0325E-10 
HH=1 ,05459- 9E-27 
K=1.380622E-16 
M0.-9.109558E-26 
P]=3. 141593 
CEVE=1.6021917E-12 
CALL FFCTS 
NHON=l 

20 IF(N~UN.E0.1> go to 21 
WHITE (6.3) 

3 FORMAT (//'CONTINUET 1»YES» 2»NO^» 

REA0<2*) NCONT 
IF (NLUNT,EQ,2) go TO 550 

21 NR0N=NRUN*1 
WRITE (6»4) 

4 FORMAT ( 'ENTER NUMBER 
6' (PAHM NO, ,NFW value 

READ(2») NCHNGEidSB 
00 23 I=1,NCHNGE 
JsISOSCP( I) 

25 PARMiJlsVALUEd) 

X2=t AhM(6» 

Xl=, ,- a2 

WRITE (6,6) 

6 FORMAT ( 'Parameters for this run arem) 

DO 2/ 1*1,7 
27 WRITE(6,) I,PARM(I) 

CALL HCSP0L(X2,THTRAN1»THTRAN2,THLI,THL2»ETRAN1,ETRAN2i 
1E0,EINF,F1RS,E2RS) 

WRITE (6, 10) E0,EINF,E1RS*E2RS,TMLIiTML2»THTRAN1,TMTRAN2# 
IETRANUETRAN2 

10 FORMAT(1X,«EPSO«*»F8.4»' EPS1NF« • ,f Bt4/IX» 'HGTE CALLEN^f 
I' CHARGE*' ,F8. 4, • COTE * " 

2' phonon temp*' ♦F8.4»( 

3'HGTE TO PHONON TEMP«tfF8.4. 

41X,*HGTE transverse CHARGC«' lF8*4f • COTE TRANSVERSE •• 
5'CHARGE ••tfe.A//) 

WHlTt(6i5) 

5 FORMAT (//• T Z CO(T» (CF»EG)/KT»/ 


I**, crsinr tr DtWAiA* .nuic, WAkLcri'f 
)T£ CAULEN CHARGE*' *F8.4/IXt 'HGTE LO' 
>• COTE LO PHONON TEMP*' tFB.4/lXt 
>tfFS«4«» COTE TO PHONON TENP>*»F8.A/ 


»PARM(n 
■ ( 2 ) 




&• 
i* 

TMIN 

TMAX*PARM( 
0TEMP*PARM(3) 
N0*PARM(4) 

ZI*1. 

NA*PARM(S) 

ZA*1, 

QEV*1.08>0,18<»X2 

0»0EV*CEVE 

UV* 0»53 


8ILITV»//> 


NLH 
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00960 

00970 

00980 

00990 

OlOOO 

01030 

01040 

01050 

01060 


£0EV«9,5 


1*X2)*1,(. 

|(1.60*6i27*X2)»l,0E5 


• 0E5 
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E2EV«3, 

VI.»(2.52*0,54»X2>«1 
VT»(1.60*0i 
M00»127,60 
Ml»200.59 
M2all2.4 

^PH]s4. •(X2*i 12. 4 *x1«200.S9«127.60>*1*660S6551/ (6.460*0.021^2) 

ALa(6.460*,021*X2)»l,0E-a 

EDISEV=PARM(7) 

EDiS=CEVE*EOISEV 
EPHI0=CEVE*E0EV 
EPHiUcEvE’-ElEV 
EPH12«CEVE*E2EV 
-\SaElRS*ES 
?Sat2R5*ES 


MlBal.66043E-24»M00«Hl/(M00*Ml» 

M2Bsl.66043E-24*MO0*M2/(MO»*H2) 


oi: 


70 


75 


90 


01670 


VA«AL*AL*AL/4, 

WL1«K«THL1/HB 
WL2aK«THL2/M8 

*NAT»?.*PMI/(M00*XI«M1*X2*M2)/1.66043E-24 

CONTINUE 

KT»K»T 

P*e^6E-8*CEVp***®’**^*®*^^^**************“***^‘*®"**®******^’ 

DELfA«E6/D 

BETA»KT/E6 

Blal./HErA 

UP«3.*Hn*H8«EG/(4,«M0«P*P» 

CA^L^SECU(yAR<2.«Pl/AL»rBC»78V#£3) 

YBVs-YBV/KT 


CONs 2.*UP*MO«KT/(M0«H01 
^ON»CON*|QRf (CONI / (2.Sp 

DEL2*ol3*Z 
GENO*. FALSE, 

X0L=FUNC(2) 

IP(XOL.LT.O.O) GO TO 75 
6EN0»,TRUE, 

XOUaXOL 


I*PI*SQRT(BETA)) 




nOENO) GO TO 90 
jO TO 85 


XOLsXOU 
Z«Z*UELZ 
XOUaFUNC(Z) 


110 


[830 


IF(XOU.Lt.O.O) 00 TO 60 
ZU>Z 

ZL«Z-UELZ 
GO TO 95 
ZLaZ 

ZUa2*OELZ 

IF (DELZ.LE,2,E-4*A0S(Zn 00 TO 110 
OELZ-OEL7/2, 

ZaZL*OELZ 

XOMsFUNCiZ) 

IF(XUM,LT,0.01 GO TO 100 
ZU»Z 

GO TO 95 

XOLaXOM 

ZLaZ 

GO TO 95 

ZaZL*OELZ* (-XOL/ (XOU-XOLl > 

XZleFUNCiZ) 

CQNa2,*uP*M0/ (HB*H8) 

C0NaC0N«50RT (CON) •2.«ES*ES«SORT (KT/BETA) /PI 
X21alO.*BI 

IF(Z.GT.el) X21«10.*Z 
XZlaAMINl (YBCtX21) 

X22al0, 

IF(Z.LT.O.O) X22al0.-Z 
X22aAMINi(YBViX22) 

CALL SiHP(FUNC3tl./BETA*X21f 0,01 t 10*51 ItKFTEZtNtlER) 

KFTE2»C0N*KFTE2 

CALL S IMP (FUNC4t0,0*X22*0, 01* 10*51 1*KFTLH2,N»1ER) 
KFTLH2sC0N*KFTLH2 
KFTHM2a(2,*UV»M0/(HB*HB) )*«(a./2,) 
KFTHH2aKFTHH2«ES*ES*SQRT(KT»»F(-l,-Z»/PI 
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01040 

olaso 

01660 

01870 

01060 

01890 

01900 

01910 

01920 

0193O 

01940 

01950 

01960 

01970 

01980 

01990 

02000 

02010 

02020 

02030 

0e:040 

02050 

O 2 O 6 O 

02070 

02080 

02090 

O 2 IOO 

02110 

02120 

02130 

02140 

02150 

02160 
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THETAlaTMLl/T 
TMETA2*THL2/T 

GINF2aHB*HB»(KFTE2*KFTUH2)2 (2, •UP*MO»EO»E1NF) 
G0^ = ^<8•Ha•(KFTE^♦KFTLH2♦KFIHH^)/■- ' "■ “ 


'(2.*UP»MO«EO«EO) 


200 


250 


-_i00 
02210 
02220 
02230 
02240 
02250 
02260 
02270 
02200 
0C290 
02300 
02310 
02320 
02330 
02340 
02350 
02360 
02370 
02380 
02390 
02400 
024>9 

oU^oc 

02430 

02440C 

02450 

02460 

02470 

02480 

02490 

02500 

02510 

02520 

02530 

02540 

02550C 

02560 

02570C 

02560 

02590 

02591 
02600 
02610 
02620 
02630 
0264 0 
02650 
02660 
02670 


A1=E2S/E1S 
Al»Al«Al*Mlfl*WU/(M2H*WL2) 

A2eEb»EPHlO/(PI«HB*ElS*ES*VL» 

A2 ca 2»A2«M1B*VA*»<L1*UP*M0/(2«*HB»PH1) 

A3sA3*4|iHB'MlB*VA*WLl*<N01«ZX*21*NA«2A*2A)/(A.«UP»M0»KT) 

A5 = ES*EG*HB«EO?S*EOlk/HB«t)p««5/HB«M^ 
l*WLl/KT/ANAT**I»X2/2,/PI/Pl 
NPlel./ «EXP(THETA1)-1,> 

NP2«1./(EXP(THETA2»-1.» 

bi=i./beta 

THLIMsahaxI (THETA lfTHETAZl 
YLIM=10.*8J*THLIM 
IFIZ.GT.BI) YL1M«10.*2 *THLIM 
ylim=amIni <YBC»YL 1M) 

00 200 I«l.3 

DO 200 J=»1.3 
N»J-1 

call SlMPlFMNtbl lYLlMf 0.001*10fSIl*MN(l»J) tLf lER) 

CONTINUE 
DO 250 I3l,3 
M»I-1 

CALL SIMP(FH,BItYLIM,0,00lf lOvSllfHLl <1) «L»IER) 

continue 

01»MNtl.i) 

a2=HN(l»l)*MN(2t2)-MN(lt2>*MN(2tl) 
Q3»MN(l,l)«(MN{2f2)*MN(3f3)-MN(2.3»«MN(3i2n 
1 ♦MN(1 .2>4(MN<2»3>«MN(3tl»-MN(2*I «MN 3t3 ) \ 

i *HN(t*3 •(MN(2»l)*MN(3t2)-HN(lf2{*MN(3tli I 
QLl»MNhtl)*MLl (2)-MLI (1)*MN(2»I) 

QI.2«HN<lf l)*(HN(2*2)*MLl (3>-MLl (2)*HN(3«2) ) 

1 ♦MN(lt2)*(MLl (2)«MN(3fi)-MN(2»i)4MLl(3)) 

2 ♦MLlh)*(MN(2il)*MN(3»2)-MN(2»2)«MNT3fih . ^ 

G 1 «ML I ( i ) *ML 1 ( 1 ) /MN ( 1 f 1 ) ♦OLl»QLl/ (Q 2 *un •QL 2 *QL 2 / ( 03 * 02 ) 

XaKT/(PI*HB»ElS> 

SlGMAaMlB«VA«WLl*X*X*Gl/(3»*PI) 

UHr6.935E6*SIGMA/NE 

EG0PTs-Bl»KT/l,602l917E-12 


£FEGKT=Z-ei 
ARlTE(6t7j TtZ.EGOPTt 


FEGKT 

12«4(2X«lPE12,4*2XtlPE12*4) 


5S0 


F0RMAT(lX,F12,4»2XtlP 

WRITE(6»9) UH 

F0RMATn3X»lPEl2.A/> 

t»t*otemp 

|F (T.l|jTmaX) j to 60 
continue 

STOP 

END 


SUBROUTINE SECULAP(KS*E1*E2»E3) 


PEAL KTtNDtMOfNA 

COMMON^/CONSTl/ PI tPtHBtKTtNOfUVvHOfOtEOtYBCfVBVtNAtZItZA 

C0»P*P*KS*KS 

C2sD»E6 

C1»-(E6*0*C0) 

CO«-2.*D*CO/3. 

call CuBE(C2tCl»C0»El»E2»E3) 

RETURN 

END 


FUNCTION FUNC(X) 


. PI,P»HBtKX*NO«UVtMOtOfEOtYBCtYBV*NAfZIiZA 

NEfNLH«NHHtN01 


PEAL NE»NLHiNHH 
PEAL KT>NOtMO*NA 
REAL NDO*NDl 
COMMON /CON5T1/ 

COMMON /NNN/ ■" 

COMMON /VAP/ ZfCON 
COMMON /C0NST2/ BETA«DEl.TA»Bl 
external FUNCI 1 FUNC 2 
NEaO.O 
X2slO.*Bl 

IF(Z.GT.BI) X2*10.*Z 
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OF<lG^^iAL PkC€ rS 

8^680 x-aminu.«cxp, of poor quality 

0^690 IF(X^,LE.B|) GO TO 100 

oItIO CAa^§IMP(^UNCl,BI,X2.0.0l.I0.Sn,NE,N. 

02720 100 NLHsO.O 
02730 X2=10, 

02740 inz.LT.O.OJ X2«ao.-Z> 

02750 X2SAMIN1 (T0V»X2) 

8^Z§0 IF Ui?,LE.O.O) GO TO 200 

Hill CALlJihP,funC2,O.0.X^.0.O1.1O 

lull 

lull 

lull 

lull g^®NE-'^‘.H-NMH-N01*Zl4NA*2A 

8l§?®- - end 

02840C •••••*••••••••••••••••*••••••••••••••••••••*## 

2?§fO- FUNCTION FViytZI 

O206OC ••••••••••*•••••••*••#••••••••••••»•••••••••• 

02870 |F(Y*Z.LE.«500«) 60 TO 10 

02900 60 TO 30 

02910 10 F9«1.0 

02920 GO TO 30 

02930 20 FV»1,/(EXP|y*2)*1,) 

02940 30 CONTINUE 

02950 RETUHN 

02960^ END 

02970C •••«••••*••••••••••••••••••••••••••••••••••••••••• 

02980 function FUNCMV) 

02990C 

"'*70 REAL LAMOA 

[O COMMON /VAR/ ZtCON 

6353? SeS'”'’'""*'*"'’*"’ 

?|8|?c ....5S2... ..................... 

03060 function FUNC2(V> 

03070C 

03080 NEAL^LAMUAV 

03090 common /VAR/ Z.CON 

2^}95 FUNC2»FV(YiZ)*LAMOAV(Y> 

03110 RETURN 

03j|o END 

03i30C 

03140 FUNCTION FVP(YiZ) 

03150C 

jF«VvZ,LE,-200.) GO TO 10 

?liJ? lUiSfs'-'-"'’'’-' ” 

?|J|? 10 

03210 GO TO 30 

03220 20 FVP>tXP(Y4Z) 

03230 FVPsFVP/(FVP4l.)/(FVP*l#) 

mi 

03260 END 

03270C ■>••••••••••••••••••••••••• 

3^2980 *•••••••«••••••••••• • 

03300 REAL LAMpA ••••• ••••••••• 

03310 COMMON /VAR/ ZtCON 

lull gUNC3=LAMOA(Y».FVPIY.-Z» 

ojsoc 

FUNCTION FUNC4(Y> 

03370c ****•******•*•••••••*••#•••••••• 

03380 real LAMOAV 

03390 COMMON /VAR/ ZtCON 

03400 FyNC4»LAMOAV<Y>*FVP(VtZl 

03410 RETUHN 

03420 END ^ 

^.^*REAL,FUNCTI0N*LAM0A(Y» ****************************** 

lit?? S£I.,DELI..ei 
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OF POOR QUALITY 

0^490 LAMDA»2t*okTA*^2LTA«0Y»BY*8r*OEL‘^4*(3.-OELTA»«By«BY 

035OO 1 ♦4,*<lc-0£UTA) •BY/3,-2, /3# 

tAM6A23^«LAMOA»SQ^f <BY« (BY-l , > • <0ELTA»BY* I . » / (X»X«J(«X»X> » /2, 

UJ5J0 RtIUWN 

03540 END 

03550C •#•*••*•**••*•••••*•*••••*••*••••••••••••••••••*•••*•••• 

03560 real function LAMOAV(Y) 

03570C ••*•••*••*••*•••*•••********•••*•••••••••••••••••••••••• 

03560 COMMON /CONST2/ BETA,DELTA»8I 

03590 BYsBETa'Y 

0360 0 LAM0AV»2.*0FLTA*0ELTA*BY*8Y*8Y-0ELTk''(3,-DELTA)»BY*BY 

03610 1 ♦4.*n,-0ELTA)*BY/3.«2./3, 

03620 X*3,**UELTA*BY/2.-1. . _ 

33630 LAMDAVx3,«laM0AV*?OPT (BY* (BY* I , ) • (0ELTA»BY-1 , I / (X*X»X«X»X) »/2i 

0J64C PF.TUHN 

1-3650 end 

03660C ••••**•*•*•*•**•*••••**'•*••••**••*••••••••••••••••••••• 

03670 function S(Y) 

03680C •*••*•••*••*•••*•••**••**••****••••*••••«’••••••••••••• 

03690 COMMON /CONST2/ BETA,DELTA»BI 

03700 BYs8ETA*r 

03710 IF(Y,F0.9I) BY«1. 

03720 S»SG)WT(BY*(8y-1,)*(DELTA*8Y*1.J/(3,»0ELTA»BY/2.*1,)> 

03730 RETURN 

03740 end 

03750C ••*•••*•*••*•*••••••••*••••••**•••••••••••••••••••••••• 

03760 real Function lkyj 

03770C 

03700 COMMON /C0NST2/ BETAtOELT A»BI 

03790 L1=S(Y)/SQRT(BETA) 

03800 y«Li*Ll*H 

03810 RETURN 

03020 END 

03eJoc #•*••*#•»*•*•••*•••*••••••*•••••••••#*••••••••••••••••• 

03840 FUNCTION H(Y, THETA) 

0385 OC *•*••••*•••*••*•*•*♦*#*•**••••*•••••••••••*••••#••••••• 

03860 H»0,0 

03870 IF(Y.LT. THETA) RETURN 

03880 H*1,0 

03090 RETURN 

03900 END ■ 

039 1 OC ••••**•*•••••••***•*••••••••••••••••••« »••••••••#•••••••• 

039*0 SUBROUTINE NABC ( Y»A»B,BP,C*CPfrU 

0393 OC •••••••••••••••••••••••••••••••*•••••*•••**•*••••**••••* 

03940 real N 

03950 COMMON /C0NST2/ BETAtOELTA, BI 

03960 QYsfablA*Y 

03970 IF(Y.EO.BI) BY«1, 

03980 DBY=UELTA*gy 

03990 0BY23srjBY*2,/3, 

04000 N»SGRT (BY* (DBY*1 .) •OGY23*2.*<BY-l,)/9,*(8y-l,)«OBY23»0BY23) 

04010 ASSORT (8 Y* (0BY*i.)*DBY2J)/N 

04020 BsS0RT(2.*(HY-1.) ) /(3,*N) 

04030 C=SQRT(0Y-1,)*DPY23/N 

0*04 0 CPsSORT ( jf3.*nflY/2,*l,)/(BV*(0BY*l,) ) ) 

Q40S0 RPsS -RT (2. ) •CR/(3.*N) 

C4(,o0 CR*0BY23«CP/N 

04070 FLs J.*(2,«D9 y»DBY*BY*0ELTA*(3.-0ELTA)*BY*BY*4»*(1,-0ELTA) 

04060 1 •aY/J.-2./3,)/(2.*(3.*OBY/2.«l,)»(3,*D0Y/2.*l#)) 

04090 RETURN 

0*100 END 

84 } ^ 0 ^ ****SUBRUUT!NE*RSL*YtRPl *RP2tRMl,RM2tsPl»SP2»SMltsM2tLBl* 

04130C •*••••**••••••••*••••••••*••••••••••••••*#•••••••••••*•••♦• 

04140 REAL LAMUA 

04J50 real lb 

04160 COMMON /CONST2/ iETA»OELTA*Bl 

04170 Common /C0NST3/ THETAl*THETA2fOINF2tO02,EPHl0tEPHn*CPHI2f 

0*160 1 VL,VTtSA2fSA3,SA4*SAS 

04190 ymT8»y-th|iai-bi 

0*290 t!J?e*y-LMETA2-0l 

U4210 sY=5(Y) 

04220 SP1*S(Y*ThET41) 

04230 SP2sS(Y*THETA2) 

04240 IF (YMIB.I-T.O.O) 60 TO 10 

04250 SM1*S(V-TheTA1) 

04260 10 IF (YM2B.(,T.0.0) 60 TO 20 

04270 SM2=S (Y-IHETA2) 

04200 20 CONTINUE 

0*290 call NA6C(Y,A.B,BP,C,CPiFL) « 

9*300 CALL NAHC(Y*TurTAl,APl*0PTf8PPl,CPl,CPPlfFLPl» 

0*310 call N»8C (Y*THE VA2,AP2fBP2#8PP2,CP2»CPP2fFLP2l 
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04J20 

0<*330 

0<»3*0 

0*350 

0*360 

0*370 

0*300 

0*390 

0**00 

0**10 

0**20 

0**30 

0***0 

0**S0 

0**60 

0**70 

0**60 

0**90 

0*500 

0*510 

0*520 

0*530 

0*5*j 

0*550 

0*560 

0*570 

0*580 

0*590 

0*600 

0*610 

C*620 

0*630 

0 * 6*0 

0*650 

0*660 

0 * 6/0 

0*680 

0*690 

0*700 

0*710 

0*720 

0*730 

0*7*0 

0*750 

0*760 

0*770 

0*780 

0*790 

0*800 

0*810 

0*820 

0*830 

0 * 8*0 

0*850 

0*860 


)488( 

0*890 

0*900 

0*910 


)*93t 

0*9*0 

0*950 

0*960 

0*970 

0*980 

0*990 

05000 

pSiS 

05030 

050*0 

05050 

05060 

05070 

05080 

05090 

05100 

05110 

05!2o 

05130 

051*0 

05150 


IF(YM)b.Lr.O.O) GO TO 30 

*0 ^onVin5£''‘^’^”^^*2#*M2.BM2.8PS2,CM2»CPM2»FLM2) 

S2=SURT(2.) 

8C=B«CPl*C«flPl 

BC»fl«CMl«C*RHl 

BC=b*CM2«C*8M2 

0a=BP"bPPi«cp*CPPi 

8C»0P»CPPl*CP*BPPi 

8X»BP«BPP1 

C*-CP»CPP1 

RlPl=2.*AiAPl*0B 

8B=|p2BPP2!cP*C?^r*®^'^^*^**®**‘^*'®^*®‘^'2.*CX.^^ 

BCeBP*CPP2«CP*8PP2 

BX=0P*BPP2 

CX«CP»CPP2 

R1P2=2,»A»AP2*BB 

IF (YMIB.LT.O.O) go to 70 

BB=BP*BPH1*CP*CPM1 

BCeBP*CPMl*CP*8PMl 

0X»BP*BPMI 

CX*CP*CPPl 

R1MI=2,*A*AM1*B0 

B8=BP*BPM2*CP*CPM2 

BC=BP*CPM2*CP*8PM2 

BX=BP*PPM2 

CX=CP*CPM2 

R1M2=2.*A*AM2*BB 

an ?2M2=3,.8 x*BX/*.*BX»BC/S2»2.«BX»CX-BC»BC/2.*CX*CX 

OW CONT INUE 
SY2=5Y*SY 
QlP=bPl*5Pl*SY2 
02 PsSP2«5P2*‘>Y2 
IMYMIB.LT.O.O) 

QlMaSMl*SHl*SY2 
90 IF (YM20 .Ct.O.O) 

32M=SM2*5M2*SY2 
100 CONTINUE 
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60 TO 90 
60 TO 100 


P1P=01P*G 

P2P*02P*6 


IFIYMIB.LT.O.O) 


no 

120 


NF2 

NF2 


NF2 
• OtO) 

NF 2 


GO TO 
GO TO 


110 

120 


o<;H=UirM*Gi 
CONTINUE 
VlPaSY'SPl 
V2P»SY*SP| 

IF'YMIBjLT.O.O) go TO 
VlM»SY*SMi 

1F(YM2b,LT, 0.0) 60 TO 
y2M=SY«SM2 
1*0 CONTINUE 

X«AL06(A8S< (P1P-2,*V1P)/<P1P*2.«VIP1 ) 1 
PVePlP«PiP-*,*VlP*VlP 


130 


130 

1*0 


.•P1P/01P»»X/*. 


V0Pl*-(2-*6lNF2*VlP/PV*X/2.»/QlP 
viPl»VlP*(P}P/PV-l,/QlP-PlP*PlP/( 01 P*Py)l*il.- 2 .»l 
V2Pi=(VlP*(I,-2.*PlP/QlP*n.»PlP/QlP)*PiP*PlP/PV> 

1 ♦P|P*(2#-3,*PlP/QlP)*X/2,>/*. 
V3Pl»-VlP*vIP*VlP/(3,*QlP»«PlP*yiP 

1 •il*-3i»PlP/(2.*01P»*PlP*PlP/<2,*PV)-PlP»PlP*PlP/(2,*01P*PVn 

2 /2,fP1P*P1P*(3.-*,«P1P/Q1P)«X/16, 

X»ALOG(APS( (P2P-2,*V2P)/(P2P*2*»V2Pn) 

PV=P2P*P2P-*,»V2P*V2P 

V0P2=-(2.»GINF2*V2P/PV*X/2.»/Q2P 

V1P2=V2P*(P?P/?v-1,/02P-P2P*P2P/(02P»PVMf( 1,-2.«P2P/02P>*X/*, 
V2P2= ( V2P* ( I . -2, *P2P/02P» 1 1 i-P2P/Q2P) •P2P»P2P/PV) 

I *P2P*(2,-3.«P2P/02P)«X/2.I/*, 

V3P2=-V2H*V2P»V2P/<3.*Q2P/ ♦P2P»V2P 

1 •< l .- 3 .* P 2 P /( 2 ,* Q 2 P )* P 2 P * P 2 P /( 2 .« PV )- P 2 P » P 2 ?* P 2 P / l 2 .* Q 2 P»PVn 

2 /2,*P2P«P2P*n,-*,»P2P/02P)»X/l6. 


146 



0b)60 

osin 

obieo 

0bl90 

0b?00 

0b210 

0b220 

0b?30 

oe>2*o 

0b250 

0b26O 

0b270 

(jbZeO 

0'52'JO 

0S300 

05310 

05320 

05330 

053*0 

05350 

05360 

05370 

05380 

05390 

05400 

05*10 

05*20 

05*30 

05**0 

0b*5Q 

05*60 


05*70 

05*8q 

05*90 


ORIGINAL PM'3£ 13 
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If IVMlH.LT.O.o: GO TO )50„,„ , 

XaAlUG ( A«S ( ( P lM-2 , • V 1 M J / ( P 1M*2 ,*V IH) ) ) 

PVaPlM*PlM-*,«VlH«VlM 

) ♦ ( 1 .-2.*P1M/01M) •X/6t 

1 ♦PlM*(2,-3,»PlM/QlM»«X/2.*/*. 

2./2.iPlM*plM*(3,-*.*PlM/OlM5«*/l6, 


150 IF (YM2B,LT,0,0» go to 160„ „ , 

X = ALOG( A3S( CP2M-2,*V2 ^1/ (P2H*2,«V2M> ) » 

PV = P2‘'«Pcm-*.»V2M*V2M , 

V0M2a-<2,*6INF?*V2M/PV*X/2,»/Q2M 
VlM2 = V2Mi (P2M/PV-1 ./Q2M-P2M»P2M/iQ2M*PVn*(l. 
V2M2= (V2M* ( 1 ,-2,«P2M/02M» ( I .-P2M/02M)*P2M«P2M/PV) 
1 *P2M»(2.-3,*P2M/02M)*X/2,»/*. 

V3M2 = -v2m«V2H*V2M/ (3,*0?M> *P2M*V2H 


2.*P2M/Q2M)«X/4* 


•(lo- 3 .»P 2 t/t 2 .*Q 2 M)*P 2 M*P 2 M/( 2 ,*PV)-P 2 MP 2 MP 2 M/( 2 .«a 2 M*PV)» 

/2,«P2M»P2rt*(3.-*.«P2M/02M>«X/l6, 


160 g^NTJNOE 




05500 

05510 

05520 

05530 

055*0 

05550 

05560 

U5570 

05580 

05590 

05600 

05610 

05e20 

05630 

056*0 

05650 

05660 

05670 

05680 

05690 

05700 

<•5710 

057?n 

Obr 30 
057*0 
05750 
05760 
05770 
05780 
05790 
05800 
05810 
05820 
05P30 
058*0 
05850 
05860 
05870 
C5880 
05890 
05900 

inn 

05930 

059*0 

05950 

05960 

05970 

05980 

05990 


=FL*FLP1* (5Pl*SPl*St2) 

SPI=KP1*3Y2 

RP1sRP 1«<P0P1*V1P1»RJPI*V2P 
SPlaSPl* (ROPf •V0P1 *RIpI«V 1P 
RP2 = fl,»FLP2* (SP2»iP2*SY2> 

SP2=RP2*5Y2 

rP2bHP2* (R0P2«V1P?*R1P2*V2P2* 

SP2=5P2» (P0P2»V0P2*RIP2*vIP2* 

IF (YMlfl.LT.O.Ol GO TO 170 
RMI»FL*FLH1« (SMl«SMl*SY2) 

SM1sWM1^SY2 

i::b4:!g§si:$isi:giai:;!!!ii:8§lii:!iiail 

GO TO |80 
RMlaO 


’}*S' 


*p; 

}*RJ 


!P2*V: 

^P2»Vi 


811 


170 RMlaO.O 

160 IfIyM^B.LT.O.Q) go to 190 
ftM2.FL*FLM2«(SM2«SM2*SY2) 


SM2«RM2*SY2 . . . 

RM2aRM2*(R0H2«VIM2*RlH2*V2M2*»<2M2*V3M2) 
SM2»5M2*(R0M2*V0M2»RIM2»VIM2»R2M2»V2R2I 
GO TO 200 
190 RM2»0,0 
SM2SO.O 

2®® PIoIa*a*a«a*B*B*0»B/*«-S2*8*8*B»C*2«*8*B*C»C 
PI 1=2,*A»A* (BP»RP*CP*CP) 

Pl2aJ,*BP«BP*BP*8P/*,*S2*0P*0P*8P*CP*CP»CP*CP*CP 
X«ALUGIAGS( (602»*.*SY2»/Gq2)) , 

XEH=AL0G(AHS( (GINF2**.*SY2)/GINF2> ) 

X1=SY2/ (G02**,*SY2> 

X15HS5Y2/ (G1NF2**.»SY2) 
pHlo=X-*.*Xl 
PHlEMO=XEH-*,*XiEH 


PHli=(G02.5Y2)*X-*,«Xl*(G02*3.*SY2»,, , 

PHlEHja(GINF2*SY2)*XEH-*,*xiEH*TGTNF2*3.»SY2) -ev»i«*y» 
PmI2»(-*.»SY2*<602*SV2) ♦ <G02*2.»5 t 2‘ • (3»«602*2»«SY2) •X/2« 


2,»Xl* <G02«2.*SY2)«(G02*2»*SY2; 1/2* 
PHIEH2« (-*.*SY?*«G!NF2*SY2»*«GINF 

•XEH/2 - - 

PHlaPIO 


H/2,-2,*xUh*]6inF 2*2,«SY2) 
....-PI0»PHl0*PIl*PMll*PI2*PHi2 
PH1 EH=p10»PHIEH0*PI1*PHIEHI*PI2«PHIEH2 


M§IN^2*I?«SyIJ 1 7l» 


•61NF2*2*«SY2) 


E I'EPHl l/FPHIO 
E2=EPHI2/EPHI0 

A1»A*A«E1*B*B/2.*E2*(C»C*B*8/2.I _ ^ , 

A2S-E1* lC*C*B«B/2.)*E2*<5.»C»C*0»B/2») 
A3s*,»t2*(C*C-9*B/2.) 

A4=EI*£1«B*B*C*C 

01»0*t>«B*B/*, 

B2C=B»S2»C 

0C2=B*C/S2 

03aB*B» (i|c*BcC«0*B26"2^»B?iH2C-BC2«BC21 

B*i!2.»fl*d*BC2*B2C 

05»'-d«B»82C*B2C 

B6«C«C-B*B/2, 

B6«(t3*B*C»C*B6*B6/3,l/2. 

FLY»AI*A1*2.*(A1*A2*A*)/3.*(2.«AI«X3**2»X2-A6)/2» 
I ♦*.»A2«A3/5,*A3*A3/3. . , . 

FTYsEl»El*(B!*2t*B2/3.*B3/2.*2»«B*/5.*B5/3.*B6) 
X3LAmuaTt) 

Xl»vL/vT 

PHI0l5»2t*A***.*B***,*2./3.*C*«*.< 


>4,*SQRT(2,)/3,*B**3*«C 
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100 


06QOO 
06010 
060^0 
06030 
060«0 
06050C ••• 
06060 

06070C ••• 
06030 
ObQ^Q 
06100 
06110 7 
06120 
06130 
06140 10 
06150 
06160 20 
06170 30 
06180 
06190 100 
06200 

Itlif ••• 

06230C 
''>240 
06250 
06260 
06270 
C6280 
06290 
06300 
06310 
06320 
06330 
06340 
06350 
06360 200 
06370 
06380 300 
06390 
06400 
06410 
06420 
06430 
06440 
06450 
06460 
06470 
06480C *•* 
06490 
06500C 
06510 
06520 
0653C 
06540 
06550 
06560 
06570 
06580 
06590 
06600 
06610 
06620 
06630 
06640 
06650 
06660 
06670 
06660 
06690 
06700 
06710 
06720 
06730 
067*0 

06 . 3 

06760 
06770 
06780 
06790 
06000 
06010 
06020 
06830 

06840 


I ♦4.»(n»C)**2,-4,/3.*A*»2.*ia«»2»4C«*2.) 

LB*SA2*X»X4SY2« (FLY*Xl*xI«FTr) ♦FL*FL*TSA3*PH1*SAMPH1EH) 
1 ♦SA5*X«X*SY2*PHI6 iS 


RETUMN 
END 


SUBROUTINE R00T(F.A#8«70L1 


Z=l.O 

EA=E (A) 

ir (EA.GT.O.Ol Z--1.0 

X=(A*B)/2,0 

Y=E{X) 

IE(Y*Z) 10*100*20 
A«X 

GO TO 30 
0sX 

IF(fl-A,LT,TOU RETURN 

GO TO 7 

RETURN 

END 


SUBROUTINE INTEG(F*R*A*B*V7 


V>0«0 

DtL=0,0001*(B-A) 

Y»X*OEL 

DEL“5t*0EL 

IE(Y.GT.8> Y«8 

CA' •. SIMP(F,x*Y*0.01*10*SP*VAt*N»iei 
WRi (6*7 X.Y*SP*VAL»N»I£ 

V»V* AL 

IF(Y.EO.B) GO TO 200 
X»Y 

GO TO 100 
OEL30.5o014(8>A) 

Y»H 

X«Y-OEL 
DEL=5,*CEL 
IFIX.LT.AJ XsA 

CALL SIMP(F*X*Y*0.01*10»SP*VAI.*N*1E) 

WR1TE(6*) X*Y«SP*VAL*N*1E 
V=V*VAL 

IFIX.EO.A) return 
Y=X 

60 TO 300 
END 

FUNCTION FMN(Y) 

REAL LR*NP1,NP2 

COMMON /YAR/ Z*C0N 

COMMON /NM/ N«M 

COMMON /CONST2/ BETA.OELTA*0I 

COMMON /CONST3/ THETA 1 *THETA2»6INF2»002fCPHlO*EPHntEPHl2» 


1 VL«VT«A2.A3*A4*A5 
COMMON /CONST4/ NP1*NP2* 


X1*X2*A1 


10 


X»Y-8I 

CALL RSL < Y*RP) *RP2*RH1 *RN2»SP1*SP2»SN1 »SH2*L8> 

XMaO.O 

XNsO.O 

IF(X.CO.O.O) GO TO 10 

XM=X*»M 

XNaX**N 

XMN=XM*XN 

xpin=(x*thetai»«*n 

XP2N3 (X*ThETA2»«»N 
XM1N*0.0 
XM2N*0.0 

IF(X-THETaI.LE.O.O) 


20 

30 


XM1N«(X-THETA1)»*N 

IF(X-ThETA2.LE*0.0> 


60 TO 20 
GO TO 30 


XM2Ns(X-ThETA2)**N 
CONTINUE 
FV=FV(Y*-Z» 

FMNxXl*Bl4XM4(FV(Y«THETAl*-Z)»(NPl*l,>»iXPlN»RPl»XN*SPl>/FY 
1 *FV U»THETAl.-Z»*NPl*MI¥»JH£TAl)*UMlN«RMi-XN«SMl)/FYl 

4 -LB*XMN 
FMN=*FMN«FVP(Y.-Z» 

RETURN 
END 


«>/FY» 


/FV 
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iVsU 

07870 

07880 

07890 


COMMON /C0NST2/ BCTA*OELTA»BI 
AmY*BI 


10 


XMaOtO 

inx.Eo.0,0) 

XMeX**M 


60 TO 10 


CONTINUE 
FM*-Fy 


VP(Yf-Z)«Ll <Y»«XM 
RETURN 
END 

SUBROUTINE HCSPOL(X*THTl*THT2tTHLltTHL2tE5TRltE5TR2tEPS0f 

IEPSItECl*EC2» 

ESTRla2,96 

ESTR2=2.35 

XMiaI.29b058E-22 

XM?a9.923396E-23 

A=<6.*60*0.21*X»*1,0E«B 

EPS I = 15.1 077-1 3,8823«X*9,88928*X»X-3.67133»X»X*X 

EP50=20.2a6-lS.1531*X«6.5909l*X*X-0.951826*X*X*X 

MT2al52.0-9.0*X 

WTlall8.0«12.0«X 


OofiboC *•*•*••••**•***•*•••••■ 

C6860 FUNCTION FMIY) 

06870C •••**•*♦•••••••••••••• 

06880 REAL LI 

06890 COMMON /VAR/ ZtCON 

06900 COMMON /NM/ N«M 

0691Q ' 

06920 
06930 
06940 
06950 
06960 
06970 
06980 
06990 
07500 
07510 
07520 

07530 

07531 

07532 

07533 
075*0 

VilU 

07570 

8JI18 

07600 
07610 
07620 
07630 
u 76*0 
076*1 
076*2 
76*3 
76*4 
76*5 
076*6 
07650 
07660 
07670 
07o80 
07690 
07700 
07710 
07720 
07730 
(i77*0 
07750 
C7760 
07770 
0776? 

0778, 


ES=*.80325E-10 
3=1.05*5919F-21 
8Ks1.380622E-16 


PI=3. 1*1593 
C=2,9979250E10 
TMTis2.*PI*C*H8/8K»HTl 
THT2=2.*PI*C«HH/HK*wt2 
0MT1=c'K*ThT1/H8 
0MT2«BK*THT2/H8 

F1=*.*PI/(A*«3./*.)*ESIR1*»2./XM1*(1,-X>»ES»*2. 

F2=4.*PI/(A**3./*J«ESTR2«*2./Xm2»X*£S«*2, 

S1=FI/0MT1**2, 

S2=F2/0MT2**2« 

THLM2al./2./EPSl*(EPSl«(THTl*»2.»THT2»*2.)*Sl»THTl»*2. 
1«S2«ThT2«*2.-SQRTHEPSI»(TmT2«#2.-TMT1**2,)*S2*THT2«*2. 
2-S1*THTI«*2.)**2.****THT1**2.«THT2«*2.*SI«S2) 1 
THLP2=l./2./EPSI*(EPSI*(THTi»*2,*THT|»*2#)*Si*TMTl»*2# 
l*S2*THt2**2.*S0RT«<EPSI«(THT2»*2,-THTl*;2,i*51»THT2*»2# 
2-Sl'*THTl*«2.)**2.**.*THTl*;2.*THT2«*2.»Sl*S2)> 
THLO=SORT< (EPS1*S1)/EPSINTMT1 
THL12«THLM2 
TML22aTHLP2 

[P (THT2.5t,THL0> go to 10 
rKLJ2sTHLP2 
THL22»THLM2 
10 THLI=S0RT(THL12) 


THL2»SQRT(THL22) 
TML1«(TML1*THL2)/2. 

/EP< 

/EP« 

RETURN 
END 

0*20C*^^*^'^SU8ROUTINE CUBE <BtC»n.El tE2»E3) 
0*21 OC 



0*?2CC 

0*230C 

0*2*CC 

0-25CC 

0*26CC 

0*270C 

0*280C 

0*29CC 

0*300C 

0*310C 

04320C 

0*330C 

0*600 

05200 

06200 

07200 

C8200 

O9200 

11200 

11300 

11*00 

U500 


CUBE 

COMMENTS i 


WRITTEN 9/19/7* 


JMP 


T«IS routine solves the general CUBIC EQUATION 

X**3 ♦ B*X«*2 ♦ C*X ♦ D ■ 0 

having only real SOLUTIONS ^E1 > E2,»„|3. 

IT IS ASSUMED THAT ALL SOLUTIONS ARE REAL* 


real E(3.' 

COMPLEX w»Xl,x2»X3,FfSO»YfZ»S 
WaCMHLX(-,5tSQRT(3,)/2.) 

Pa-(0*B/3.) *C 

Q=(2.*B*B*B/27.»-8C*B/3.)*0 

F=!0*Q/*.)*(P*P*P/27.) 

SO = CSORT <F) 

s=sc-a/ 2 . 

Y = 0,0 

IFiCABSlS) .NC.0.0) V«CEXP(CL0G(S)/3.) 
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11900 

S=-5U-Q/2, 

12000 

Z=»0,0 

12100 

If (CABS(S) .NE.0,0) : 

1h200 

XlsY.Z-B/3, 

15200 

X2=W*Y*W*W*Z-B/3« 

16200 

X3=W4W4V*W*Z-B/3. 

16300 

E(U*REAL(Xll 

16400 

E(2)»REAL(X2) 

16500 

E (3 “REAL(X31 

16600 

00 100 l»l,2 

16700 

J»l*l 

16800 

DO 100 K«J,3 

16900 

IF(E(1>,GE.E(K)I go 

17000 

»*E ( I) 

17100 

E(l: =‘F(K) 

1 7200 

E(K)*X 

17300 100 

continue 

17400 

E1=E(1) 

17500 

E2«6(|) 

17600 

E3aEl3) 

16200 

RETURN 

19200 

END 


Z«CCXP(CL0e(S)/3.) 
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TO 100 


_ ?CEI INTRPLtMO306B 
SUBRrjTINE fPCTS 


00050SrRFt 

Oooeossouf 

00100 SUBF . - . , 

00120 COB> jN ETATA012*!) fFTM12(24l),FTl2(2*n fFT32l2*l)»FT52l2*l»t 

00125f-FT-'2(2'l) 

00130 r',!_L nRFN(4, fFRTBLS.M0326B »#0tn 
DO 10 I'l»24l 

00150 10 ETATF3(I)s<I-4n/10. 

>lt241)f (FT12(I)fI*lt2«l)tlFT32mtI«lt2«l) 

FT72n> tI«lt2Al) 




00160 READ(4. > (FTM12( I ) t 
00170^*lFT52i:>»I=l*24l)f 
001 75 CALL CI.OSE (4) 

00180 HETURN 
OOiso END 

00200 FUNCTION F(K2»ETA) 

00210 DIMENSION C0N(4»5) 

00220 COMMON ET AT A8 ( 241 > ,FTM1 2 (241 ) tFTl2 (241 ) 
00230 data C0N/1.7724539fl.4l42136f-.4il2335l.- 
00 24 0*1, 88622693 » 2, 028427 1 ♦1.233700S»l.0654n9t 
00250M.329 34 04.,5,6568542»6,l685O27t-l.7T56866» 
00260*3.323351,11 .3137081 14. 393i73tl2.429806* 
^711*111.86825/ 


1 *7756866* 


00270*11.63 1728, 22,6274 17 ,25, 907 
00280 Ik=(K2,3)/2 
00290 X=ETA 

00300 IF<X,GE,-4,> GOTO 10 

00301 EXPXsO.O 

00302 IF (X.LE.-200.) GO TO 5 

00310 EXPX=EXP(X) 

00311 5 CONTINUE 

0U320 F=C0N(l,IK)*EXPX*(l,-EXPX/C0N(2*IKn 
00330 RETURN 

00340 10 IF(X.GT.20.) goto 20 

00350 CALL INTRPL (X ,£T AT A0,FT12 ( 1 )• 17.241 *F, INO) 

00360 IF(1N0,£Q.0) RETURN 
00370 WRITF(6,8001) K2,X 

00380 80C1 FORMAT ('INTERPOLATION ERROR IN F ( « * 12* */2. • *012*5« * ) • ) 

00390 CALL EXIT 
00400 20 XM2*X*«(-2) 

00410 XKPl»IK-,5 

00420 F>X**XKPl/XKPl*(l.«CONOtIK)*XH24CON(4tIK)*XH2*»2) 

00430 RETURN 
00440 END 

OOIOOC^LAGPANGIAN INTERPOLATION ROUTINE 
00110 subroutine INTRPL(X.TA8X,TABy*NCtN*Y»lNO) 

‘0120C x=indepenoent variable 

TABX=X array 
TA8Y=Y array 


8oJ30C 

00140C 

00150C 

00160C 

00180C 

00190C 

00200C 

00210C 


NC - CONTROL WORD 
TENS POSITION 
1 Np^ 

0 YES 

UNITS POSIT 
OF values 


N»NO, 
Y = 


--^lOl 

DE^EWENT‘’VA9IA8LE*BEfyRNE0,p^CALLING PRO^^ 

O0P20C iNO=ERROR INDICATOR RETURNED TO CaCCTnG PROGPAI 
00230 DIMENSION TaBX (I) ,TA0Y (1) 

00240 ;'C=IABS(NC> 

00250 U<KC.GE.10) GOTO 65 
00260 KEY-1 
00270 GOTO n 
00280 55 KEY*! 

00290 KC=KC-10 
00300 70 IOX=KC 


ON 


EXTRAPOLATION OUTSIDE RANGE OF X IN TABLE 
DEGREE OF INTERPOLATION (1«7) 


r 
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ERROR 


0U310 1ND=0 
C0320 XA=X 

00330 CALL OISSER<XA,TABXtNtIOXtNPXltKEYtlNO) 

003^0 CALL Lft'jUAN(XA,TABX(NPXl) .TA8Y(NPXU »I0 a*a. I. 

00350 HETURfJ 
00360 Eno 

00370 SUriROUTINE OISSER (XAtTABtNX. lOtNPXtKEY, !N0) 

0O330C nlbSER LOCATES XA IN TAB AND RETURNS ilTM 
00390C NPx=lOwEW LIMIT FOR SUMMATION.' IN0«l IF E 
OOaOO OI'-'ErJSION TaB(1) 

00410 IF(NX.ST.ID«1) GOTO 70 
00420 ID=NX-1 

00430 70 IF(XA-TAB(1)| 71.73.72 
00440 71 IhOsiND*! 

004S0 IF (KEY.EO.O) GOTO 7A 
00460 73 XA=rAaii) 

00470 74 NPX»1 
00..80 RETURN 

OC490 72 IF(XA-TAB(Nxn 1.77.76 
OOSCO 76 INC=1ND.1 
OOSIO IFIKEY.EQ.01 GOTO 78 
00S20 77 XA5TAH(NX) 

OOSJO 70 NPX=NX-IO 

OOS40 return 

00550C IF DI53ER TAKES THIS PATH XA IS WITHIN RAN6C OF TAB 

C0560 1 NPT=1D*1 

00570 IF (N. .FO.NPT) GOTO 7A 

00500 NF'h = l;Pi/2 

00590 NPU = ,NPT-NPB 

00600 10 NL0w=l*NP8 

00610 NUDPrNX-NPU 

00620 IF (Nx.LE.LO) GOTO 18 

00630 NXX=;jX/2*l 

00640 IF (XA-TaB (NXX) ) 12.17.13 

00650 13 NXX=NXX.NX/4 

00660 IF(aa.GE,TAB(NXX) I GOTO 17 

00670 NL0K=NXX-NX/4 

00600 GOTO 10 

00690 12 NXX=NXX*NX/4 

007G0 iF(XA.LT.TA8(NXX)l GOTO 18 

00710 17 NL0W»NXX 

00720 18 00 19 lUNtOW.NUPP 

00730 NLOC«II 

00740 IFITaBUI) .GE.XA) GOTO 20 

00750 19 continue 

00760 NPX=NUPP-NPB«1 

00770 return 

00700 20 NL=NLOC-NP0 

00790 NU=Nl*I0 

CoSOO 00 25 JJaNL.NU 

OLAlO NDISsJJ 

U0P^Y 1FiTab(JJ).EQ.TA8(JJ41)) GOTO 30 
00030 25 CONflNUE 
OC0«»O NPX^NL 
OL05O RETURN 

00060 30 !F(TaB(NDIS) .LT.XA) GOTO AO 

OOA/o NPx=riCl5-I0 

00000 RETURN 

0O890 40 NPXsNDlS*! 

00900 RETURN 

00910 end 

0U920 subroutine lagranua.x.y.n.sumi 

00930C LAGDaN performs the summation 

00940 dimension Xdl.Ydl 

00950 SUMaO.O 

00960 DO 3 Isl.N 

00970 PROO«Y(I> 

00900 DO 2 Jsl.N 
00990 AsX(I>-X(J) 

01000 IF(A.EQ.O.O) GOTO 2 
OloiO B»(XA-X(J) )/A 
01020 PR0D*PH0D*8 
01030 2 CONTINUE 
01040 3 SUM>SUM4PR00 


50 RETURN 
60 END 


SUBROUTINE SIMPlF.A.B.OELtlMAX.SIl.SfN. ICR) 


OOIOOSFIXED 

oonoc * “ 

00120 
0OJ30C ' 

00140C SIMPSON INTEGRATION ROUTINE. 

00170 N30 

00180 0A=B-A 

00190 IF (BA120.19.20 
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